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PREFACE 


The  electrical  behavior  of  a  semiconductor  is  critical  to 
its  usefulness  in  a  broad  spectrum  of  applications  ranging  from 
infrared  sensitive  detectors  and  microwave  amplifiers  to  high 
density  integrated  circuits.  This  behavior  is  controlled  by 
extremely  low  concentrations  of  impurities  -■nd  defects  in  the 
semiconductors  and  by  its  temperature.  The  only  method  available 
for  quantitatively  measuring  small  impurity  concentrations  and 
determining  their  identities  is  based  upon  the  Hall  effect.  The 
electrical  transport  properties,  namely  resistivity.  Hall  mobility, 
and  carrier  concentration  are  calculated  as  functions  of  temper¬ 
ature  by  combining  resistivity  and  Hall  effect  measurements.  The 
temperature  dependent  carrier  concentration  data  is  analyzed  to 
obtain  the  concentrations  and  identities  of  the  impurities.  Con¬ 
centrations  of  less  than  5xlOi;icm-3  (1  part  in  1011)  are  measured 
routinely  in  silicon. 

Interpretation  of  the  Hall  data  is  limited  by  our  incomplete 
knowledge  of  the  Hall  scattering  factor  (or  r-factor)  which  is  the 
ratio  of  the  Hall  mobility  to  the  conductivity  mobility,  r=y,,/u  . 
The  data  analysis  calls  for  the  conductivity  mobility  while  the 
Hall  mobility  is  provided  by  the  experiment.  Uncertainties  in 
dopant  concentration  as  great  as  50%  can  result  when  little  or  no 
information  is  available  on  the  precise  values  and  the  temperature 
dependence  of  the  r-factor.  Unfortunately,  no  straightforward 
experimental  method  is  available  to  measure  r(T)  over  all  temper¬ 
atures  and  impurity  concentrations  of  interest,  due  in  the  case 
of  silicon  to  material  and  facilities  limitations,  mainly  a 
sufficiently  large  magnetic  field. 

Theory  of  charge  carrier  scattering  can  in  principle  be  used 
to  calculate  y(T)  over  all  physically  realistic  temperature  and 
impurity  ranges.  By  comparing  theoretical  functions  with  experi¬ 
mental  results  taken  over  a  sufficiently  broad  region,  parameter 
selection  can  be  verified  and  the  theoretical  functions  can  be 
extended  with  confidence  to  temperature  and  impurity  ranges 
unreachable  experimentally. 


Two  primary  developments  prompted  undertaking  a  new  and 
more  thorough  solution  to  the  r(T)  problem,  namely  the  increasing 
technological  importance  of  more  precise  silicon  analyses  and  the 
availability  of  ultra-high  purity  silicon  that  made  new  experi¬ 
mental  measurements  feasible.  A  coordinated  experimental  and 
theoretical  program  was  planned  and  implemented  to  solve  the  low 
temperature  region  in  p-type  silicon  where  charge  carrier  trans¬ 
port  is  limited  by  acoustic  phonon  scattering.  New  results  are 
presented  for  silicon  and  germanium  and  impurity  analysis  has 
been  significantly  improved  as  demonstrated  for  the  silicon  case. 

It  is  planned  to  extend  the  theoretical  calculations  to  higher 
temperatures  through  the  inclusion  of  optical  phonon  scattering 
which  will  build  directly  upon  the  acoustic  phonon  work.  It  is 
hoped  to  eventually  complete  the  low  temperature  region  for  doped 
semiconductors  by  solving  the  ionized  impurity  scattering  problem. 

The  theoretical  work  was  performed  under  Air  Force  Contract 
Nos.  F33615-78-C-5064  and  F33615-81-C-5095  at  the  Materials 
Laboratory,  Wright-Patterson  Air  Force  Base,  Ohio.  We  are 
indebted  to  Mr.  John  A.  Detrio  of  the  University  of  Dayton  Research 
Institute  for  his  long-standing  support  of  this  work.  w.  C.  Mitchel 
gratefully  acknowledges  support  by  a  National  Research  Council 
Fellowship.  We  would  also  like  to  acknowledge  Drs.  R.  S.  Allgaier, 
P.  J.  Price,  M.  Tiersten,  J.  C.  Hensel,  S.  Zukotynski,  C.  Elbaum, 

J.  E.  Lang,  and  C.  Searle,  for  many  helpful  discussions.  A  special 
note  of  thanks  goes  to  Miss  Niki  Maxwell  for  a  fine  and  meticulous 
job  in  typing  the  manuscript. 
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SECTION  I 
INTRODUCTION 

The  Hall  effect  has  long  been  one  of  the  principal  means 

of  characterizing  semiconductors.  The  hole  concentration  can 

be  determined  from  the  relation  p=r/eR  where  R  is  the  Hall 

n  H 

coefficient  and  r  is  a  dimensionless  parameter  of  order  one  called 
the  Hall  factor  which  depends  on  both  the  band  structure  and  the 
scattering  mechanisms  of  the  material  being  studied  as  well  as  on 
the  magnetic  field.  Measurement  and  curve  fitting  of  the  carrier 
concentration  as  a  function  of  temperature  yields  information  on 
the  identities  and  concentrations  of  dopants  and  impurities  in 
the  material.  In  addition,  the  carrier  concentration  can  be  used 
with  the  conductivity,  a,  to  determine  the  mobility  from  the 
definition  u  =a/pe=aR  /r.  Because  the  Hall  factor  is  not  known 

C  H 

precisely,  it  is  often  set  equal  to  one  for  convenience.  Unfor¬ 
tunately,  this  assumption  is  strictly  valid  only  for  degenerate 
semiconductors  or  when  the  magnetic  field  is  high  enough  that  the 
magnetic  motion  dominates  the  scattering.  This  occurs  when 
w  t>>1.  Here  u)c=eH/m*  is  the  cyclotron  frequency  and  x=am*/pe2 
is  the  scattering  time.  H  is  the  magnetic  flux  density  and  m*  is 
an  effective  mass.  For  many  practical  situations  the  condition, 
w  t>>1,  is  unattainable.  For  this  reason  a  distinction  is  made 
between  two  mobility  definitions.  The  conductivity  mobility  is 
defined  as  u  =a/pe  and  the  Hall  mobility  is  defined  as  u„=oR„. 

The  ratio  of  the  two  mobilities,  y„/y  is  therefore  the  Hall 

ri  C 

factor.  The  conductivity  mobility  is  the  one  most  often  cal¬ 
culated  by  theorists  and  the  low  field  Hall  mobility  is  the  one 
most  often  experimentally  determined.  This  unfortunate  situation 
has  developed  over  the  years  because,  while  the  Hall  mobility  is 
the  easiest  to  measure,  the  conductivity  is  the  easiest  to  cal¬ 
culate  since  it  contains  only  one  component  of  the  magnetocon¬ 
ductivity  tensor,  o  ,  but  the  Hall  mobility  contains  two,  a 

XX  XX 

and  a  .  A  similar  situation  exists  when  one  attempts  to  use 
Hall  effect  data  to  determine  the  dopant  or  impurity  concentra¬ 
tion,  but  here  no  formal  distinction  is  made  between  the  actual 


1 


concentrations  and  those  determined  by  assuming  a  simple  form 
for  r.  Traditionally  the  form  of  r  used  to  analyze  the  data  is 
mentioned  in  the  report  so  there  is  rarely  any  real  confusion, 
but,  unless  the  Hall  effect  measurements  were  made  in  the  high 
field  limit  or  the  correct  r  used  in  the  data  analysis,  the 
resultant  dopant  concentration  can  be  off  by  as  much  as  a  factor 
or  two  when  compared  with  the  results  of  other  techniques  such  as 
C-V*1*  or  neutron  activation  analysis (2*. 


We  have  measured  the  magnetic  field  dependence  of  the  Hall 
coefficient  as  a  function  of  temperature  of  ultra-pure,  p-type 
silicon.  We  have  determined  both  the  high  and  low  field  limits, 
wct>>1,  and  <<1  respectively.  We  have  determined  the  Hall  factor 
from  the  relation  r=R^ (0) /R^ (») ,  and  we  have  also  determined  the 
low  field  Hall  and  the  conductivity  mobilities,  all  in  the  temper¬ 
ature  range  20-50K. 


Theoretical  modeling  in  support  of  the  experimental  effort 

was  motivated  by  a  general  lack  of  theoretical  calculations  for 

p-type  silicon  carried  from  first-principles  without  adjustment 

of  parameters  to  fit  the  data.  Acoustic  phonon-limited  electronic 

(3—25)  (26—28) 

transport  in  diamond  and  zincblende  type  structures 

has  often  been  modeled  within  the  framework  of  the  deformation 
potential  theory.  Within  the  theory  the  strength  of  the  carrier- 
phonon  interaction  is  governed  by  an  experimentally  determined 
set  of  deformation  potential  parameters.  These  parameters  are 
obtained  from  experiments  which  measure  shifts  and  distortions 
of  the  bands  upon  application  of  selected  strains  to  the  crys- 
tal  .  However,  in  some  versions  of  the  scattering  theory 

it  has  been  found  necessary  to  regard  these  interaction  parameters 
as  adjustable  in  order  to  obtain  a  quantitative  agreement  with 
the  transport  data.  This  is  particularly  the  case  in  the  trans¬ 
port  calculations  for  p-type  materials ( 1 °" 14 ' 17 ' 21 ~24 *  where 
either  effective  deformation  potentials  or  the  relaxation  times 
were  fitted  to  the  transport  data. 


The  need  for  adjustment  of  parameters  arises  as  a  result 
of  approximations  introduced  in  the  theoretical  treatment  of  the 


2 


p 


full  deformation  potential  scattering  theory.  It  is  a  common 
simplification  to  use  a  relaxation  time  approximation  for  the 
solution  of  the  Boltzmann  equation.  One  of  the  most  elaborate 
examples  of  such  a  model  is  pi-j/ided  in  the  work  of  Bir, 

Normantas,  and  Picus^6).  This  model  employs  separate  intra- 
and  interband  relaxation  rates  which  depend  on  averaged  deforma¬ 
tion  potential  constants  whose  magnitudes  are  varied  in  the 
calculation.  Other  authors  have  considered  separate  relaxation 
times  for  each  band,  with  a  common  effective  deformation  potential, 
but  without  any  interband  coupling ( 10 ' 14 ' 16 ' 22 ' 25> . 


Along  with  the  relaxation  time  approximation,  the  degree 
to  which  the  full  nonspherical-nonparabolic  nature  of  the  valence 
bands  has  been  taken  into  account,  at  various  points  in  the  cal¬ 
culation,  has  varied  from  one  theoretical  treatment  to  another. 

(  6 ) 

The  simplest  approach,  as  for  example  that  taken  by  Bir  et  al . 

(  14  ) 

and  Costato  et  al .  ,  is  to  consider  just  the  two  top  valence 

bands  -  the  heavy  hole  band  and  the  light  hole  band.  Both  bands 

are  treated  as  spherical-parabolic.  A  slightly  more  complex 

model  treats  the  light  hole  band  on  the  average  as  nonparabolic. 

(16) 

This  model  was  set  by  Ottaviani  et  al.  who  also  considered  a 

single  nonspherical  but  parabolic  band  model  along  with  a  Monte 

(16) 

Carlo  technique  in  order  to  study  the  effects  of  anisotropy 

(22) 

Asche  and  von  Borzeszkowski  ,  on  the  other  hand,  studied  the 
effect  of  band  nonparabolicity  on  mobilities.  In  order  to  account 
for  band  nonparabolicities  as  accurately  as  possible  they  calcu¬ 
lated  the  hole  density  of  states  factors  in  the  relaxation  rates 
by  using  a  spherically  averaged  Kane ( 40 5  secular  equation,  which 
results  in  isotropic  bands.  Following  Asche  and  von  Borzeszkowski, 
Nakagawa  and  Zukotynski ^ 10 '  performed  a  calculation  in  the  same 
spirit  but  did  not  spherically  average  the  bands.  And,  more 

(25) 

recently,  Takeda,  Sakui  and  Sakata  repeated  the  work  on 

Nakagawa  and  Zukotynski  but  in  addition  considered  the  spin-orbit 
split-off  band's  contribution  to  the  mobility  .  In  view  of 
the  number  of  approximations  introduced  in  previous  treatments 
of  the  problem  it  appears  worthwhile  to  approach  the  problem  again 
from  a  first  principles  point  of  view,  without  the  band  shape 
or  the  relaxation  time  approximations. 
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Except  for  the  treatments  of  Tiersten^  and  Lawaetz^4* 
for  Ge,  a  complete  implementation  of  the  deformation  potential 
theory  has  not  been  pursued  from  first  principles.  These  two 
authors  did  not  use  the  relaxation  time  approximation  but  instead 
each  solved  the  Boltzmann  equation  directly.  Tiersten  used  the 
deformation  potential  theory  of  Whitf ield ^ 4 1  *  and  the  vector  mean 
free  path  theory  of  Price ^42'  to  solve  the  Boltzmann  equation. 
Lawaetz,  on  the  other  hand,  used  the  transition  rates  of  Bir  and 
Picus^5^  and  solved  the  Boltzmann  equation  by  expanding  the  dis¬ 
tribution  function  in  cubic  harmonics.  Both  methods  were  found 
to  agree  with  one  another  to  within  2%  for  the  conductivity  mobil¬ 
ity  of  Ge  and  were  within  a  few  percent  of  the  experimental 
mobility.  In  addition,  Lawaetz  calculated  the  Hall  mobility  and 
other  galvanomagnetic  coefficients.  Here  he  found  that  his  results 
could  not  be  reconciled  with  the  independently  measured  set  of 
deformation  potential  parameters.  Partly,  this  discrepancy  can 
be  attributed  to  the  incomplete  knowledge  of  the  parameters,  or 
to  the  quality  of  data  used  for  comparison  with  the  theory.  In 
our  work,  however,  we  will  show  that  the  results  for  the  Hall 
mobility,  at  least,  can  be  improved  by  treating  the  light  hole 
band  as  nonparabolic,  which  was  done  by  neither  Tiersten  nor 
Lawaetz . 

No  comparable  treatment  to  that  of  Tiersten  or  Lawaetz  has 

been  tried  for  p-type  silicon.  Whereas  germanium  has  a  large 

0.30  eV  spin-orbit  splitting,  in  silicon  the  corresponding 

splitting  is  only  0.044  eV.  Therefore,  its  bands  are  decidedly 

more  nonparabolic  than  is  the  case  for  germanium.  Even  in  the 

case  of  germanium  we  will  show  that  light  holes  dominate  the  Hall 

mobility,  and  the  light  hole  band's  nonparabolicity ,  however 

small,  does  influence  the  results.  To  this  end  we  have  calculated 

hole  —  acoustic-phonon  transition  rates  for  Si  and  Ge  using 

experimentally  determined  deformation  potential  parameters.  We 

(4344) 

pointed  out  in  two  previous  papers  '  ,  that  inclusion  of  the 

full  nonparabolic  nature  of  the  Si  valence  band  system  results  in 
a  strikingly  energy  dependent  set  of  transition  probabilities. 
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Parabolic  bands,  though,  lead  to  energy  independent  transition 
probabilities . 

The  goal  of  the  present  work  is  to  fully  implement  the 
deformation  potential  theory  for  the  calculation  of  conductivity 
and  Hall  mobilities  in  Si  and  Ge,  with  emphasis  on  Si,  in  the 
acoustic  phonon  limited  regime.  By  using  the  eigenfunctions 
obtained  from  the  numerical  solution  of  the  full  6x6  k*p  Hamil¬ 
tonian  we  make  no  approximations  concerning  the  nonspherical- 
nonparabolic  nature  of  the  Si  valence  band  structure.  Our  trans¬ 
ition  probabilities  are  energy  dependent  and  depend  on  both  the 
incident  and  scattered  hole  wave  vectors  as  well.  We  solve  the 

Boltzmann  equation  directly,  by  a  variant  of  the  method  used  by 
(  4 ) 

Lawaetz  ,  for  the  distribution  function.  Then,  analogously  to 
Lawaetz,  we  obtain  the  conductivity  and  Hall  mobilities. 

In  Section  II  we  formulate  the  transport  theory  to  be  used 
in  the  present  work.  Our  formalism  differs  in  a  few  details  from 
that  of  Lawaetz  who,  in  addition,  has  not  provided  the  full 
derivation.  Given  the  complexity  of  the  formalism,  we  include 
in  Section  III  details  of  the  computational  procedure  and  ways 
to  simplify  its  various  steps.  Section  IV  is  used  as  a  test  of 
the  computational  scheme.  Here  we  use  the  parameters  of  Tiersten 
and  Lawaetz  for  Ge  and  exhibit  an  overall  agreement  between  our 
respective  results.  In  the  case  of  the  Hall  mobility  we  show  an 
overall  improvement  in  our  results  as  compared  with  experiment. 
Section  V  is  devoted  to  the  presentation  of  experimental  and 
theoretical  results  for  Si.  Conclusions  are  presented  in 
Section  VI. 


SECTION  II 
THEORY 


Our  task  is  to  solve  the  Boltzmann  equation  for  the  hole 
distribution  functions  for  p-type  Si  when  both  the  electric  field 
F  and  magnetic  field  B  are  present.  With  the  distribution  func¬ 
tions  we  are  to  calculate  both  the  conductivity  mobility  y  ,  the 
Hall  mobility  y  ,  and  their  ratio,  y„/y  =r  —  the  so-called 

ii  rl  C 

Hall  r-factor. 


The  Boltzmann  equation  is 


(2) 


_|-F  •  B  •  (v+e 

*h  k  -ft  \  k  N 


(iox^) 


fN(k>  =  SopfN(k) 


(1) 


where  f.T(k)  is  the  steadv-state  distribution  function  for  valence 
N  ^  ^ 

band  N  of  energy  EN(k)  at  wave  vector  k  in  the  Brillouin  zone  of 
Si.  SQp  is  the  scattering  operator  given  by 


SoPVk)  = 


V 


(277) 


X  /’dk'[fM(k')-f 

M=1  J 


N(k)]WNM()c,k'),  (2) 


V  is  the  volume  of  the  crystal,  and  the  sum  over  M  is  carried  over 

the  three  top  valence  bands  of  Si.  The  transition  probabilities 

for  hole  scattering  by  acoustic  phonons  between  bands  N  at  k  and 

(44) 

M  at  k'  were  expressed  in  a  previous  paper  as 


Here  kfi  is  the  Boltzmann  constant,  T  is  the  absolute  temperature 
and  (k)  is  a  normalized  cubic  harmonic  in  direction  k  with 
angular  index  L  transforming  like  the  yth  row  of  the  ith  irre¬ 
ducible  representation  of  0^  with  dimension  l The  coefficients 

A  were  fitted  to  the  transition  rates  in  a  previous  paper  on 

(44) 

angular  variations  of  the  rates  .  These  rates  were  calculated 
in  Ref.  43  using  the  deformation  potential  theory  of  Tiersten(3). 
In  that  paper,  the  three  top  valence  bands  of  silicon  were  used 


taking  into  account  both  their  nonparabolicity  and  anisotropy. 

Similar  rate  calculations  were  performed  by  Tiersten  and  Lawaetz 

in  the  parabolic  band  regime  for  Ge^3'4^.  In  the  rather  good 

elastic  scattering  approximation  of  holes  by  acoustic-phonons, 

the  expansion  coefficients  A  are  only  functions  of  one  of  the 
(44) 

holes*  energies  .  The  delta  function  in  (3)  reflects  this 
energy  conservation  condition. 


as 


(4) 


The  Boltzmann  equation,  1,  can  be  written  more  concisely 


lf-fop  +  S-Sop,fN(it>  -  SopV*>  • 


(4) 


where  from  (1) 


<%V*)x?e) 


(5) 

(6) 


To  first  order  in  both  the  electric  and  magnetic  fields,  the 

(  4 ) 

steady-state  distribution  function  can  be  written  as 

fN(£)  =  +  F^0(k)  +  FB^fk)  ,  (7) 

N 

where  reflects  the  deviation  of  the  distribution  function 

from  its  thermal  equilibrium  value 


h.  -(E  (k)-u)/k  T 

f°(k)  =  e  =  f£(EN(k))  ,  (8) 

where  y  is  the  chemical  potential,  upon  application  of  the  electric 
N 

field.  similarly,  is  the  added  deviation  from  f^(k)  when 

additionally  a  low  magnetic  field  B  is  turned  on.  Inserting  Eq. 

(7)  for  fN(Jc)  into  Eq.  (4),  and  equating  the  same  powers  of  F 
and  §  on  both  sides  of  (4) ,  one  obtains 

*  Wio**1  • 

and 

*  SoAl(h  ' 

where  F  and  B  are  unit  vectors  in  the  directions  of  F  and  B, 
respectively.  Equations  (9a)  and  (9b)  reflect  the  ascending 
hierarchy  of  equations  needed  to  evaluate  high  order  FnBm 


(9b) 


perturbations  of  the  distribution  function.  Equations  for  general 
order  n,m  in  the  series  have  been  exhibited  by  Lawaetz (4).  For 
calculations  of  conductivity  and  Hall  mobilities  we  need  go  no 
further  than  the  first  two  equations  9a  and  9b.  This  is  consistent 
with  typically  low  magnitudes  of  the  electric  and  macmetic  fields 
used  in  the  Hall  experiments.  Incidentally,  the  hot-electron 
effects,  in  the  non-Ohmic  regime,  necessitate  higher  order  F  powers 
in  the  fN(k)  expansion,  Eq.  (7). 

We  solve  the  Boltzmann  equation  using  the  method  developed 
by  Lawaetz.  For  this  we  expand  and  in  a  cubic  harmonic 

series 


■  fN(*>  S  • 

-  fs<*>  S  x5;,<v*>>k!;?(£> 


(10a) 

(10b) 


with  expansion  coefficients  <p  and  x»  respectively.  Note  that  the 
most  rapid  energy  variations  of  and  4^  have  been  factored  out 

in 

f°(k)  =  f°(EN(k))  ,  (11) 

leaving  relatively  smooth  functions  <J>  and  x* 


1.  SOLUTION  OF  THE  BOLTZMANN  EQUATION 

Our  method  for  solving  Eqs.  (9)  is  sufficiently  different 

from  that  of  Lawaetz  to  merit  its  own  development  in  this  sub- 
.  (  2 ) 

section.  Whereas  Lawaetz  multiplies  both  sides  of  Eqs.  (9) 
on  the  left  by  a  cubic  harmonic  (k)  and  integrates  over  the 

angle  k,  we  in  addition  multiply  these  equations  by  6(EN(k)-E) 
and  integrate  over  the  whole  Brillouin  zone.  The  field  term 
in  Eq.  (9a)  becomes  then 

F*  fdk  K^(k)6(EN(k)-E)^f°  (St)  .  (12) 
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In  order  to  carry  out  the  Briilouin  zone  integral  in  Eq.  (12), 
we  define  for  each  band,  N  =  S,  L,  and  H  (spin-orbit,  light,  and 
heavy  valence  bands,  respectively),  functions  N,  where 

n  2  m  a 

kN  E  YN(EN(k),k)  ,  (13) 

and  m  is  the  f ree-electron  mass.  Most  conveniently 
o 


2m 


=  1  Ry  ( a .  u . ) 


(14) 


where  Ry  =  Rydberg  =  13.6058  eV,  and  1  a.u.  =  atomic  unit  of 


length  =  1  Bohr  radius  =  0.5292  x  10 


-8 


cm. 


As  we  will  show,  the 


atomic  units  are  a  more  natural  set  of  units  to  carry  out  the  cal- 

(43  44 ) 

culation  than  the  units  employed  in  two  previous  publications  ' 


Eq.  13  suggests  the  change  of  the  integration  element  to 


7  i/2mo\TT  *  dYM(E„(ic)  ,k> 

k2dk  =  -g)2Y^(EN(k)  ,k)~  - dE„(k) 


-h' 


dEN(k) 


N 


(15) 


Altogether,  with  Eqs.  (8)  and  (15),  Eq.  (12)  assumes  the  form 

'r'  (16) 


ef°(E) 

2fiik  T  F‘^NX'  (E) 

d 


where 


/ 2m  \|-  r  A  a  y 

5ft.<E>  =  (^f)  /dk  K^,(k)Y2(E,k) 


dYN(E,k) 

dE 


(V^EN(k))E,  (17) 


and  the  energy  gradient  is  evaluated  on  the  constant  energy  E 

_  2 

surface.  Since  y  is  in  Rydbergs,  g  is  in  (a.u.) 


The  collision  term,  on  the  right  hand  side  of  Eq.  (9a) ,  can 
be  similarly  simplified.  But  first  we  note  that  the  expansion 

coefficients  A  in  Eq.  (3)  for  W„„  were  given  in  Refs.  41  and  42  in 

2  3  NM  . 

units  of  (eV)  cm  /erg.  This  prefactor  can  also  be  written  in  atomic 

units  as  7.9537  x  1011  Ry(a.u.)3.  In  the  old  units  the  A's  were  or. 

the  order  of  10  11  so  that  in  the  new  units  the  A's  are  on  the  order 

of  unity.  We  will  keep  the  atomic  units  until  comparison  is  made 

with  the  mobility  data.  Multiplying  Eq.  (9a)  by  (k) 6 (E^ (k) -E) 

-41 

and  integrating  over  k  and  k',  with  given  by  Eq.  (10a),  SQp 


1 


i 

I 

i 


V 
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defined  by  Eq.  (2),  and  V?NM  by  Eq.  (3),  and  the  change  of 
variables  Eq.  (15),  we  obtain 


R.H.S . 


1  kBT 

1  - B -  f  (E) 

4  (4ir)*fi  N 


£ 

kMnl ' 


sjk 

NvX",MnX 


(18) 


Here,  we  define  the  scattering  matrix  S  as 


(F)  ■= 

t>NvX",MnX' 


E 

iULL' 


<jvX»|r(J|il,L>AjtiML,<iUL'|rH|knX-> 
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NM 


z 


(19) 


1 J  <jvX»,knX' irN|iUL>A^LfPL,<iyL'  | Tp 1 1 > , 


The  bra-kets  denote  angular  Brillouin  zone  integrals  between 
normalized  cubic  harmonics 

<k  j  j v X >  =  K^v (k)  ,  (20) 

and  the  local  operator  TN  is  defined  as 

-  -  -  /2m  \|-  \  -  dyM ( E , k ) 

rN(E,k)  =  <k|rN(E)lk>  =  ^j2Y2(Efk)_N__ -  .  (21) 

It  is  very  important  to  the  rest  of  the  analysis  that 

TN(E,Rk)  =  TN(E,k)  ,  (22) 

A 

for  all  rotations  ReO^,  the  cubic  point  group,  so  that  TN(E,k) 
transforms  like  the  identity  representation.  The  property  (22) 
follows  trivially  from  the  defining  Eq.  (21)  for  TN  and  (13) 
for  yn/  where  one  has 

EN(Rk)  =  EN(k)  (23) 

from  the  invariance  of  the  electron  Hamiltonian  with  respect  to 

space-group  operations  of  the  crystal.  It  is  easy  to  show  that 

3  - 1 

the  collision  matrix  S,  Eq.  (19),  has  units  of  (Ry  a.u.  ) 

The  advantage  of  the  form  Eq.  (19)  over  the  one  used  by  Lawaetz 
is  that  our  matrix  S  is  symmetric  in  indices  (NjvX,  MknX'). 


Since  TN  belongs  to  the  identity  representation  we  may 
employ  orthogonality  of  the  cubic  harmonics  with  respect  to  the 
representation  and  partner  indices;  i.e.,  in  Eq.  (19) 


j=k  and  v=n 


10 


so  that 


.jk 


,  (E)  -  S-l  „  m,,(E)(S..6  ,  = 


NvX"  ,  Mn  X '  ‘  NX",  MX'  jk  vX 

-L 


LL 


'  <jvX"| rN| jvL>A^L>ML,<jvL' I rM| jvX 
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NM 


<jvX",  jvX'  |  r  | 

p=i 


aL>ASL,PL'<CtL' 


(24) 


We  use  oi  in  Eq.  (24)  to  denote  the  one  dimensional  identity 
representation . 

Combining  results  of  Eqs.  (16),  (18),  and  (24)  the  following 

system  of  linear  equations  results: 


-2e 


sjv 

NX" ,MX ' 


<p 


jv 

MX' 


(25) 


Eq.  (25)  is  to  be  solved  for  the  expansion  coefficients  <J> .  As  can 
be  seen,  the  system  of  linear  equations  factors  out  into  separate 
problems  for  different  irreducible  representation  j  and  partner 
indices  v.  If  we  go  one  step  further  we  can  show  that  solutions  $ 
have  a  simple  dependence  on  the  partner  index  v.  From  the  matrix 
element  theorems  we  can  write  the  integrals  in  Eq.  (24)  as 
follows  * 45 * 


i  . 


<jvX"|  rN|  jvx> 


l 

i . 
i 


t 

u=l 


j  U  X  ‘ 


N 


j  U  x  > 


(26) 


for  all  l<v<*j,  where  l ^  is  the  dimensionality  of  the  jth  irre¬ 
ducible  representation  of  0^,  since  TN  transform  like  the  identity 
representation.  Similarly,  the  "Gaunt-like"  integral  of  three 
cubic  harmonics  in  Eq.  (24) 


< jvX", jvX' | Tn | aL> 


x,  . 

=  ^  ^  <jyx",j 
j  u=i 


yX' |rlaL> 


where  TN|aL>  transforms  like  the  identity  representation  since 
<k | aL>  is  a  cubic  harmonic  for  the  identity  representation.  Eqs. 
(26)  and  (27)  show  then  that  the  index  v  on  S  in  Eqs.  (24)  and 
(25)  is  redundant. 


Examining  the  g  vector,  Eq.  (17),  shows  that  j  must  belong 
to  the  vector  representation,  r^,  since  ^E^  (k)  transforms  l  ike 
r*15  (or  6  in  the  notation  of  Von  der  Lage  and  Bethe  -  VdLB46)  . 
For  all  other  j  the  L.H.S.  of  Eq.  (25)  is  zero  so  that  4>=0  for 
j^<5.  Rows  of  the  6  representation,  v,  can  be  labelled  by  the 
Cartesian  coordinates  x,  y,  and  z.  The  matrix  element  theorem 
expressed  by  Eq.  (26)  can  be  used  for  g  as  well,  with  the  result 
that  for  fixed  v 


^NA"  (E)  = 


.  r  ~  -  /3E  (k)\ 

5“  £?Jdk  KX”(k)VE'k> 


["<GNX" 


(E)  /  2)  e. 


where  e^  is  a  unit  vector  in  the  vth  direction.  Using  the  above 
developments,  Eq.  (25)  becomes 


,(E)*iw(E) 


NX  ",  MX '  ' ^MX 


which  has  the  solution 


bNX-  "  e(kBl)  (F'ev)  (SNX'  .MX"  <E))  GNX" 

=  e(^)2(?.;v,  [eiA,m]  , 


where  the  directional  dependence  is  contained  in  F*e^.  If  we  def 
define 


VE-k> 5  ?  •(#)' 


KXX(^>«x 

K5xZ(k)ez 


then  to  first  order  in  the  electric  field  Eq.  (7)  can  be  written 
in  the  familiar  form 

fN(ic)  =  f»(lt)  [l+?.$N(E,k)  1  .  (31) 

Clearly,  the  units  of  and  4>  are  those  of  the  (electric  field)  1 
while  those  of  6  are  (Ry  a.u.). 


Eq.  (29)  provides  the  solution  for  that  part  of  the  distri¬ 
bution  function  which  is  perturbed  by  the  application  of  the 
electric  field  only.  The  second  part  of  the  problem  is  the 
solution  of  Eq.  (9b)  for  the  part  of  the  distribution  function 
which  is  perturbed  by  the  simultaneous  application  of  the  mag¬ 
netic  field.  As  for  the  R.H.S.,  we  use  the  expansion  (10b)  for 

•  A 

4^,  multiply  the  R.H.S.  by  (k)  6  (EN  (ic) -E)  and  integrate  over  k. 
This  part  of  the  calculation  is  identical  to  the  calculation  above 
with  x  replacing  tp,  so  that 


R.H.S.  =  =- 


kBT 


(4tt) 


fN(f) 


Z 

MX' 


’NX"  ,MX ' 


(E)X 


3^ 

MX  ‘ 


(E) 


(32) 


where  x  are  the  expansion  coefficients  for  41!:  ..  in  Eq.  (10b)  .  On 
the  L.H.S.,  Eq.  (9b),  4-^g  is  the  solution  obtained  earlier  in 
Eq.  (30),  and  the  is  given  by  Eq.  (6).  ®op'  due  to  the 
presence  of  the  cross-product,  commutes  with  any  function  of 
energy,  such  as  f^(E  (k)),  matrices  S  and  G,  and  expansion 
coefficients  4.  Multiplying  the  field  term,  Eq.  (9b),  by 


*{?.  (k)  6  (En(£)-E) 


and  integrating  over  k,  we  obtain 


L.H.S.  =  — ~f 
2fi2 

where  the  magnetic  field  matrix  is  given  by 


jUU,jv  ,F . 

NX'", NX"  [t} 


-I 


dk  K^(k)TN(E,k) 


[%EN(J)X?j]EKi“(k) 


(33) 


(34) 


evaluated  on  the  constant  energy  E  surface.  The  sum  over  yX"' 

in  Eq.  (33)  comes  about  because  of  the  cubic  harmonic  expansion 

N  -1 

of  4^g.  Clearly,  index  i=6  and  u=x,y,  or  z.  H  is  in  (a.u.) 

As  we  shall  show  later,  Eq.  (52b)  for  the  Hall  mobility  cal- 

N 

culation  we  need  only  the  part  of  4^  which  transforms  like  the 
vector  representation  ^^5=6,  so  that  j=5  as  well.  We  shall  not 
calculate  the  part  of  4^  transforming  like  the  rePresenta“ 

tion  since  it  is  irrelevant  for  our  goals.  We  note  that 
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■Vb’ 


transforms  like  a  pseudovector  according  to  the  F'  =<5'  representa¬ 
tion.  The  magnetic  field  matrix,  Eq.  (34),  can  be  shown  to  be 
nonvanishing  for  the  choice  i=j=S  since  the  direct  product 


ri5®rI5  ’  ri5+r25  +  ri  +  ri: 


contains  in  its  decomposition.  The  pseudovector  operator  h 

transforms  under  coordinate  rotations  according  to 


P  (R)  h  P  { R" 
P 


rm  (r) h 
qp  q 


where  P(R)  are  the  scalar  rotation  operators  which  operate  on 
functions,  as  opposed  to  rotation  matrices  R  which  operate  on 
vectors.  rm(R)  is  the  matrix  for  the  pseudovector  F '  irreducible 
representation  of  0h 

With  these  preliminaries,  we  can  show  that 
”NV'\NA"  (E)  =  <6vA"[rNh(E)  |6yA'">  (38) 


EE, 

ReOh  p  y 

z  z 


e  <P(R)<5vA"|  r  (E)P(R)h  (E)P(R  l)  |  P  (R)  |  <5u  '  A'"  >  . 
p  1  N  p 


“i A  ,  Vv'v(R,rqp<R)rS'P(R)<6V,X"|rN(E)hq(E)|6U'X'''>- 

ReO^  p , q , v  ,  y 

Summing  over  the  48  rotations  of  0^,  we  obtain  the  following  result 
rt6y,6v  ,  r ,  (39) 


■tZ 

p 

-iZ 


\  7  i  e  € 

6  p  vpy 


=-  /Lj  e  € 

6  p  vpy 


,  ,  ,  vpy 

v',P',y'  ^ 


[WNA"',NA"  (f). 


,  <  6  v '  A  "  |  F  (E)  h  ,  (E)  ]  5  y  '  A"'  > 


N\  —  7  t 

P 


where  €  is  the  Levi-Cevita  symbol, 
vpy 
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Therefore,  the  L.H.S.  in  Eq.  (33)  becomes 


L’H'S*  12fi2  fN(E)  ^NV"(E)  (ep*B)€vpyHNA'",NA"  (E)  '  (40) 

and  substituting  Eq.  (29)  for  <J>  gives 


L. 


■  ■  ■  M  ^5 


z 


xm  eNA"'(E)  (F'ey)  ^^(ep‘B)evpyHNx»'fN>»  • 

P  (41) 


Summing  over  u  and  p,  we  finally  obtain 


S*kBT/  12  " 


0NA"' ^ E ^NA'" ,  NA"  ^E ^  ^FxB^  "ev  ^  '  (42) 


In  combination  with  the  R.H.S.,  Eq.  (32),  the  linear  system  of 
equations  has  the  solution 


5  v 
(PA 


/  4t 

\3  4-rre2 

[\kBT, 

*  J 

N^"(SpA,NX") 

3 


A  A 


em-(E>HN\"’,N  x",E» 


K  I(FxB,.evl  ,  |(A'' 


B 


(FxB-eJ  ^X(E) 


Defining 


(43) 


XN(E,k) 


E 


4 it  V  47re2 


A  L\kBTy 


•fi 


4<E> 


f  A  A 

K,  lk)e 

Kjy(k)e 
A  y 

Sz  "  " 

KA  (k)ez 


(44) 


the  distribution  function  to  first  order  in  the  magnetic  field 
becomes 


fN(k)  =  f°7(k)  [l+F*fM(E,k) +(FxB)  *X„(E,k)  ] 


N 


N 


'N 


(45) 


2. 


TRANSPORT  COEFFICIENTS 


The  average  thermal  velocity  in  a  material  is  given  by 


<V>  = 


IE 

k  N 


vN(k)fN(k) 


EE 


(46) 


k  N  fN(k) 


where 


v*>  --sW*1 


In  Appendix  A  we  derive,  for  completeness  sake,  the  relationship 
between  the  current  density  J  and  the  fields  to  be^47^ 

J  =  pe<V>  =  o  F+a  u„  (FxB)  +  ...  ,  (47) 

O  On 

where  p  is  the  carrier  density,  o  is  the  conductivity 

°o  =  peuc  '  (48) 

Uc  is  the  conductivity  mobility,  and  uH  is  the  Hall  mobility.  Hall 
experiments  measure  uH,  although  the  data  analysis  for  carrier  den¬ 
sities  requires  y  .  These  are  related  by  the  so-called  r-factor 

r  =  UH/PC  •  (49) 

We  shall  develop  the  necessary  theoretical  framework  for  evaluating 
Uu  and  U  • 

n  C 

Since  the  distribution  function  contains  terms  proportional 
to  F  and  FB,  we  can  write  Eq .  (46)  in  the  following  form 

<V>  =  <V1q>F+<V11>FB  ,  (50) 

which  through  Eqs.  (47)  and  (48)  provide  us  with  the  means  of  cal¬ 
culating  mobilities.  These  equations  are 

U  c  =  F*<V10>  (51a) 

ryc2  =  (FxB)  -<Vi;L>/ |  FxB  |  2  .  (51b) 
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Performing  the  fields  of  expansion,  Eq.  (50),  we  obtain 


uc  ■ 


kN 


fVk) 


E-E  (( FxB )  • 


ru 


2  _  kN  pX 


((FxB)  *VN  (ic>)  f=(i:)X^(EN(E)jKju 


(k) 


(52a; 


( 52b  ] 


FxB 


kN 


fN<k) 


It  is  clear  now  that  due  to  the  vectorial  nature  of  VN(k)  only 
the  5  parts  of  fN  contribute  to  the  mobilities.  Denominators  in 
Eqs.  (46),  (52a)  and  (52b)  are  the  total  number  of  carriers  in 

the  crystal. 

The  integrals  occurring  in  Eqs.  (52a  and  b)  can  be  written 
in  terms  of  integrals  evaluated  earlier.  The  forms  we  adopt  are 


uc  -•*  F 


E  f 

12,  NyX  Jo 


dE  e  ^NA  E  WE) 


i.E  / 

N  JQ 


dE  e  ^E<a  |  F  |  a  > 
o  1  N  1  o 


where  6  =  (kgT)-1,  and 


2  _  FxB 
r  u  —  — * — * - 

I FxB I  2 


I  / 

NpX  Jo 


) 

-BE  5y  ->-6u 
dE  e  *NA{E)gNA(E) 


4tt 


E  / 

N  Jo 


dE  e  ^E<a  |  r._  |  oc  > 
o  1  N 1  o 


(53a) 


(53b) 


Since  (F*e^)  and  g0uve^,  Eqs.  (29)  and  (28),  respectively. 


►<5vl 


Uc  ^  (P*eu)  =  1 F  i  2  =  1 

so  that  pc  is  independent  of  the  magnitude  and  direction  of  the 
electric  field,  as  it  should  for  ohmic  mobility  of  cubic  materials. 
Similarly,  from  Eq.  (43) 
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so  that 


X5uMFxB) 


ruc  = 


^  or—  ^ 

(FxB)  «eyJ  |(FxB)  *ey 
- - - * - * - =  1  , 


FxB 


and  uH  is  independent  of  the  directions  and  magnitudes  of  F  and 
as  long  as  they  are  not  parallel. 


The  symmetry  analysis  presented  in  this  section,  with  some 
modifications,  can  be  carried  over  to  noncubic  materials.  A  main 
difference  will  arise  due  to  the  fact  that  the  three  Cartesian 
coordinates  will  no  longer  be  equivalent,  and  as  a  result  mobilities 
will  be  anisotropic. 


In  the  next  section  we  will  show  how  to  implement  the  above 
developed  formalism. 
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SECTION  III 

COMPUTATIONAL  APPROACH 


In  order  to  implement  the  theory  outlined  in  Section  II  we 
need  to  evaluate  a  number  of  Brillouin  zone  integrals.  In  all 
cases  the  weighting  factor  is 


defines  yKT(E,k).  The  y„  for  all  three  bands  are  obtained  from 
N  N 

the  solution  of  the  cubic  equation  as  given  by  Madarasz,  Lang, 

,  tI  (48) 

and  Hemenger 

H3y3+H2y2+H1y+Ho  =  0  ,  (55) 

where  the  coefficients  H.  are  functions  of  energy  and  the  angle 

/N  1 

k.  Band  parameters  A,  B,  and  C  for  Si  are  taken  to  be  dimension¬ 
less  so  that  y's  come  out  in  Rydbergs.  It  is  a  simple  matter, 
given  a  general  solution  of  Eq.  (55),  to  take  its  energy  deriva¬ 
tive  to  form  the  TN  for  all  three  top  valence  bands. 

One  type  of  angular  integrals  which  may  be  performed  with 
relative  ease  are  those  contained  in  the  definition  of  the 
scattering  matrix.  The  simplest  of  these  is 


■  / 


dk  K^(k) 
Li 


U2 )  ve- 


dYN(E,k) 


All  the  functions  in  the  integrand  transform  like  the  identity 
representation,  so  that  the  integral  can  be  performed  over  one 
irreducible  wedge  of  the  Brillouin  zone  and  the  result  multiplied 
by  48  -  the  number  of  such  wedges  in  the  Brillouin  zone.  The 
numerical  technique  we  employ  is  that  of  12  x  12  point  Gaussian 
integration;  i.e.,  12  points  for  the  azimuthal  angle  and  12  points 
for  the  polar  angle  in  Eq.  (56) (10’48^  _  As  a  check,  one  may  note 
that  in  the  parabolic  limit 
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E(k) 


(57a) 


at  E  =  0.00450  eV  to  be  8 . 0 3x10 - 3 ( Ry -1a . u. _3 )  for  the  H-band  and 
1.16x10  (Ry  a.u.  )  for  the  L  band.  By  comparison,  using  the 
nonparabolic-nonsphe'rical  solutions  for  y  from  Eq.  (55)  ,  the 
computer  calculates  at  this  very  low  energy  8 . 08xl0-3 (Ry_1a . u. ~3 ) 
and  1.17x10  3(Ry  1a.u._3)  for  the  H  and  L  bands,  respectively. 

The  small  difference  between  these  numbers  gives  us  confidence  in 
our  ability  to  perform  these  integrals,  and  points  to  the  fact  that 
the  parabolic  approximation  at  very  lowest  energies  is  a  fair 
approximation  to  the  true  bands.  All  angular  integrals  calculated 
here  were  also  calculated  using  a  7x7  Gaussian  mesh  with  little 
difference  from  results  of  the  larger  12x12  mesh. 


Another  type  of  integrals  encountered  in  the  evaluation  of 
the  scattering  matrix  S,  Eq.  (26)  ,  involves  cubic  harmonics  for 
the  6  representation. 

3 

<6vA"| rN| 5vA>  =  j  ^  <6yX"| rN| 5pX>  (59) 

y=i 

which  in  fact  are  independent  of  the  partner  index  v.  in  order 
to  reduce  Eq.  (59)  to  an  integral  over  one  irreducible  wedge  we 
derive  the  following  identity: 


f 


l . 

3 

r-  2  <juL|r  |juL'>  = 
j  M=1 


rIS  /  kl 

j  u-l  RcO. 


rW(Rk)  rM(E,Rk)KT3V  (Rk) 

L  N  L 


h  W 

where  W  denotes  one  irreducible  wedge,  and 


rN(E,Rk)  =  TN(E, k) 

(Rk)  =  P(R_1)K^U(k)  = 


E  -j(R)uXJ<k> 


so  that  from  the  greater  orthogonality  theorem 


(45) 


l  . 
3 


<jULirNi:ui->  - f 

3  y=i  3  u=i  •'w 


dk  K3P(k)  r(E,k)K3? (k)  .  (60) 

L  N  L 


For  the  identity  representation  2.^=1,  which  explains  the  factor 
of  48  in  connection  with  Eq.  (56).  For  the  representation 

the  premultiplier  is  16  since  2.^  =  3. 

The  third  type  of  angular  integrals  encountered  in  the 
calculations  is  that  for  g  in  the  field  term  of  the  Boltzmann 
equation,  Eq.  (17) , 


5jK  -  <J'>»lr,l<«S»B(S))E> 


(61) 


Since  the  energy  gradient  transforms  accordingly  to  the  vector 
representation,  j=6 ,  Eq.  (61)  can  be  reduced  to  the  form  of 
Eq.  (60)  with  £.=3.  In  order  to  obtain  the  energy  gradient  we 
differentiate  the  secular  equation,  Eq.  (55),  with  respect  to 
Cartesian  components  of  k.  After  some  tedious  algebra,  one  can 
show  that 


(sirjj.  -  V  K‘E- 


k)+z3n2(E,k)+i-n3(E,k)J 


(62a) 
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where  expansion  coefficients  IL  ,  l<i<6,  are  too  lengthy  to  repro¬ 
duce  on  these  pages.  The  important  point  is  that  the  K' s  were 


constructed  so  that  they  transform  according  to  the  identity  rep- 

9EN  3  1 

resentation,  therefore,  rrr— ,  z  ,  ~  has  three  factors  transforming 

d  K  2 

Z 

like  the  vector  representation,  r^^=5.  As  a  check  on  this  part  of 
our  computational  procedure  we  note  that  in  the  parabolic  limit 


"z  'E 


z(*>  (#b2  fr  d*  =  <6zllrNl  (w~ 


>  = 


z  E 


2mo\ ( m*  \  8tt_p  /  3 

\*2)\mo)  3  *47T 


(63) 


l~T~ 

y  4tt 


which  for  the  H  band  at  E=0. 00450  eV  is  7.85x10  4a.u.~2  and 
-4  -2 

2.17x10  a.u.  for  the  L  band.  Computer  generated  numbers  for 

-4  -4  -2 

these  two  cases  are  7.45x10  and  2.14x10  a.u.  ,  respectively. 

The  fourth  type  of  integrals  is  the  most  involved.  It 
arises  in  the  evaluation  of  the  magnetic  field  term  vector,  H, 
Eq.  (34).  The  H  matrix  is  given  by 

hn“"'!nv’,e)  =  K^(k)rN(£,k)  [?ev£>*»fl  Kx“<1'1 


/ 


dk  K^(k)TN(E,k)hN(E,k)  K^,(k) 


(64) 


where  hN  operator  was  defined  in  Eq.  (35).  The  angular  integral 
can  be  reduced  to  a  one  wedge  integral  from  the  known  transforma¬ 
tion  properties  of  the  factors  in  the  integrand 


HNA"',NX" 


■E..E  E  .  rj. 

q 

x  L 


_  _  .  ,  (R)rm  (R)r5,  (r) 

q  ,  v  'v  qp  y  y 

q  M  ReOh  v  , p, y 


(65: 


w  dk  kJ?/  (k)  rN(E,k)hq(E,k)KjU,'  (k) 
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Summing  over  the  48  rotations  in  0^  yields 


u ,  6  v 

MNA",NA' 


e  € 

p  puv 


(p'u  '\ 


’puv 


dk  K 


6  v 

A" 


(k) rN(E,k) 


hp,  (E,k)K^/  (k)  j  .  (66) 

In  practice,  we  choose  the  direction  of  the  electric  field  along 
the  x-direction  so  that  u=l,  the  magnetic  field  in  the  v-direction 
so  that  p=2,  then  we  only  require  the  v=3  part  of  x  in  order  to 
calculate  the  Hall  mobility. 


We  still  have  to  calculate  the  bracketed  term  in  Eq.  (66). 

We  find  it  most  effective  to  write  out  the  magnetic  field  oper¬ 
ator  h  in  the  form  of  Eq.  (35). 

hN(E,k)  =  (V£En(£)x$£)e  , 

with  the  energy  gradient  expanded  in  the  form  of  Eq.  (62).  We 
then  form  all  the  necessary  gradients  of  the  cubic  harmonics  in 
order  to  evaluate  the  bracketed  term  in  Eq.  (66).  The  resulting 
expressions  are  tedious  to  derive,  are  lengthy,  and  inelegant. 
Therefore,  we  will  not  repeat  them.  Another  method  for  evaluating 
the  bracketed  term  in  Eq.  (66)  is  suggested  in  Appendix  B. 

It  should  be  said  at  this  stage,  for  benefit  of  other  workers 
wishing  to  perform  similar  calculations,  that  in  spite  of  our 
emphasis  on  employing  the  cubic  harmonics  the  entire  calculation 
can  be  performed  without  them.  The  lengthy  expressions  in  eval¬ 
uating  the  derivatives  of  cubic  harmonics  are  partly  the  result  of 
their  orthonormality.  One  can  form  symmetry  adapted  cubic  functions 
which  have  correct  transformation  properties  yet  need  not  be  ortho¬ 
gonal  to  one  another.  As  an  example,  for  the  5  representation  we 
could  have  used  simple  functions  like  x,y,z;  x3,y3,z3;  x5,y5,z5; 
etc. ,  rather  than  the  orthogonal  set 


x 


x 


3 


x 


3 


3 


5* 


10 

T 
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This  would  simplify  the  gradient  expressions  considerably.  That 
the  use  of  cubic  harmonics  is  not  essential  is  due  to  the  form  of 
the  angular  integrals  we  are  evaluating.  All  of  them  contain  a 
weighting  factor  between  two  cubic  harmonics.  Therefore,  we  never 
have  a  chance  to  use  the  orthonormality  of  cubic  harmonics  with 
respect  to  their  angular  momentum  indices.  The  same  was  true  in 
our  discrete  least-squares  fitting  of  the  transition 
probabilities ( 44  * . 

As  always,  we  tested  this  part  of  the  calculation  in  the 
parabolic  limit.  For  example,  in  the  parabolic-isotropic  limit 

/2m  \I  II 

y*dk6lz(k)r(k,E)hy(k,E)6lx(k)  =  2-fi^fj2  '  <67) 

which  for  the  L-band  at  E=0.0045  eV  yields  0.0145  a.u.  1  and  for 
the  H  band  0.0277  a.u.-1.  The  corresponding  results  for  the  actual 
silicon  valence  bands  are  0.0144  and  0.0268  a.u.  1,  respectively. 

With  all  the  angular  integrals  calculated  we  proceed  as 
follows.  Using  the  expansion  coefficients  A  for  the  transition 
rates  we  form  the  S  matrix  in  Eq.  (24).  The  X ',  X" indices  ranged 
through  angular  momenta  1,3,5,  and  5'  (for  the  6  '  harmonic). 

The  angular  momenta  for  the  a  cubic  harmonics  were  L=0,4,6,  and  8. 
The  S  matrix  inversion  was  performed  using  the  Gauss-Jordan 
method ^ 49 *  *.  As  noted  earlier,  our  form  of  the  scattering  matrix 
S  is  symmetric  in  the  combined  NX  indices.  The  symmetry  of  S  also 
serves  as  a  check  of  our  computational  procedure.  In  the  next  step 
we  multiply  the  $  vector  by  the  S  1  matrix  and  calculate  <J>  in 
Eq.  (29)  .  Similarly,  we  calculate  the  x  coefficients  from  <J>  and 
H  on  L.H.S.,  Eq.  (33),  and  the  S  matrix,  as  indicated  in  Eq .  (43). 

The  last  sensitive  procedure  is  the  evaluation  of  the  energy 
integrals  (53a  and  b)  containing  the  Boltzmann  factor.  Here,  we 
use*  the  Gauss-Laguarre  integration  method  with  25-points  along  the 
energy  scale We  perform  a  change  of  variables  so  that  the 
coordinate  of  the  25th  point  occurs  at 

E25  =  -(kBT> lnO.OOi  (68) 
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We  also  tried  E^ ^=-  (kgT) InO . 01  with  little  change  in  the  results. 
This,  of  course,  required  us  to  interpolate  the  integrand  at  the 
adjustable  Gauss -Laguerre  mesh.  The  interpolation  was  done  via 
the  Lagrange  interpolation  using  16  first  principles  points  as 
input (51). 


All  the  energy  dependent  quantities  were  evaluated  at  two 
separate  energy  meshes.  The  "rough"  energy  mesh  had  the  energy 
interval  of  0.0045  eV  and  16  energies  from  0.0045  eV  to  0.072  eV. 
The  same  quantities  were  also  calculated  on  a  finer  mesh  near 
E=0  eV  with  the  mesh  size  of  0.0135/17  eV;  again,  16  energies  with 
that  mesh  size  were  used.  Using  the  finer  mesh  we  fitted  the 
energy  dependence  of  the  numerators  and  denominators  in  Eqs. 

(53a  and  b)  in  a  convenient  form  for  the  interpolation  purposes. 

It  was  our  observation  that  for  energies  less  than  the  first  mesh 
point  the  Lagrange  interpolation  was  ineffective,  so  the  fitting 
was  indicated.  From  presence  of  the  Boltzmann  factor  it  is  clear 
that  the  integrals  derive  much  of  their  magnitudes  at  these  low 
energies . 


The  energy  dependence  for  the  integrands  can  be  obtained 

from  a  simple  analysis.  At  low  energies  the  expansion  coefficients 

A  for  WMN  in  Eq.  (3)  are  essentially  energy  independent,  as  if  to 

a  first  approximation  the  bands  were  parabolic.  This  makes  the 

2 

scattering  matrix  proportional  to  f  ,  so  that 


S^F 


(69) 


From  (63) ,  we  also  see  that  g^E,  so  that 


4>^S  '*'g^^Evconst . 


(70) 


which  implies  that  the  numerator  in  Eq.  (53a)  behaves  like  E  at 
low  energies.  Even  though  this  is  a  polynomial  form,  the 
necessity  for  fitting  the  energy  dependence  of  the  numerator 
arises  from  the  general  weakness  of  the  Lagrange  interpolation 
near  the  endpoints  of  the  supplied  mesh. 
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(71) 


Next,  we  observe  that,  Eq.  43, 

x^s"1$« 

Eq.  (67)  shows  that  H^E^ ,  so  that 

,  (72) 

and 

xg^  •  (7  3) 

u 

Therefore,  the  numerator  in  Eq.  (53a)  goes  as  E  and  has  to  be 
fitted  to  that  form  for  low  energies.  Similarly,  from  Eq.  (58), 
we  see  that  the  denominators  in  Eqs.  (53a  and  b)  go  as  E^.  The 
square  root  dependence  explains  why  the  Lagrange  interpolation  is 
ineffective  for  energies  less  than  the  first  mesh  point.  The 
Lagrange  interpolation  is  a  polynomial  interpolation  scheme  and 
would  have  a  difficult  time  fitting  a  function  like  E'5  near  end¬ 
points  of  the  energy  mesh. 

In  Sections  IV  and  V  we  shall  display  both  the  numerators 
and  denominators  of  Eqs.  (53a  and  b)  for  both  germanium  and 
silicon.  We  consider  the  calculation  for  germanium  using  Tiersten's 
and  Lawaetz's  parameters  as  a  test  of  the  overall  theory  developed 
in  Section  II  and  the  computational  procedure  developed  in  this 
section.  The  test  is  provided  in  the  next  section. 
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SECTION  IV 


TEST  OF  THE  THEORY  AND  COMPUTATIONAL 
PROCEDURE  -  CALCULATION  FOR  GERMANIUM 


The  only  two  calculations  similar  in  spirit  to  the  present 
one  were  performed  in  the  60' s  by  Tiersten (3)  and  Lawaetz *4)  for 
germanium.  In  order  to  assure  ourselves  of  the  validicy  of  our 
theory  and  computational  procedure,  we  performed  the  calculation 
for  germanium  using  the  same  set  of  input  parameters  for  Ge  as 
did  Tiersten  and  Lawaetz.  In  practice,  this  entails  a  change  of 
about  10  lines  in  our  code  for  silicon. 

Table  1  provides  the  parameters  which  serve  as  an  experi¬ 
mental  input  into  this  part  of  the  calculation.  We  have  to 
emphasize  that  there  is  a  large  spread  in  experimentally  measured 
values  of  the  deformation  potential  parameters  for  Ge.  From 
tabulations  of  various  measured  values  for  a,  b,  and  d  for  Ge 
it  is  apparent  that  there  is  a  large  uncertainty  in  their  mag¬ 
nitudes  and,  ac  times,  in  their  signs  as  well ( 4 ' 28 ' 35 ' 52 *  . 

Since  this  part  of  the  calculation  serves  mainly  as  a  check  of 
our  computational  procedure  we  adopt  the  Tiersten’ s  set  of 
parameters.  This  way  we  assure  ourselves  of  a  one-to-one  com- 
parison  with  the  earlier  works  in  the  field.  Moreover,  Tiersten 
and  Lawaetz  did  motivate  their  selection  of  the  deformation 
potential  parameters  from  physical  considerations.  We  then 
expect  that  these  values  are  reasonable  representations  of  the 
experimental  state  of  affairs  in  measurements  for  Ge. 


TABLE  1 

EXPERIMENTAL  INPUT  PARAMETERS  USED  IN 
THE  MOBILITY  CALCULATION  FOR  GERMANIUM 


Valence  Band 
Parameters3 
(dimensionless ) 


Deformation 
Potentials3  (eV) 


Phonon  Parameters3 
oc2  (dyne/cm2) 


A  -13.27 
B  -  8.63 
C  12.4 


a 

2.0 

s=transverse 

5.75x10 

b 

-2.1 

s=longitudinal 

1.53x10 

d 

-7.0 

3  Ref.  3,  p=mass  density. 


cg=phonon  speed  for  branch  s 


With  the  set  of  parameters  in  Table  1  we  have  calculated 

the  transition  probabilities  for  Ge  in  the  manner  described  in 

our  first  paper  on  the  subject*43*.  These  were  then  fitted  to  a 

double  cubic  harmonic  series  using  the  method  treated  in  detail 

for  Si  in  our  second  paper*44'.  The  expansion  coefficients  A, 

Eq.  (3),  that  we  calculated  for  Ge  are  in  closest  agreement 

with  those  reported  by  Tiersten  for  transitions  involving  the 

(  4  ) 

more  parabolic  H  band.  In  turn,  Lawaetz  reports  that  the 
coefficients  he  calculated  also  agree  with  those  of  Tiersten*3*. 

Before  presenting  the  final  mobility  results  for  Ge,  it 
would  be  worthwhile  to  exhibit  some  of  the  intermediate  results 
for  eventual  comparison  with  Si.  These  intermediate  results 
relate  to  the  solution  of  the  Boltzmann  equation  and  have  not 
been  displayed  by  either  Tiersten  or  Lawaetz.  But  first,  we 
should  keep  in  mind  a  few  differences  between  our  respective 
calculations .  We  use  the  three  top  valence  bands  by  solving  the 
cubic  equation,  Eq .  (55),  for  the  valence  band  dispersions.  As 

a  result,  our  bands  are  slightly  nonparabolic,  whereas  Tiersten 
and  Lawaetz  use  only  the  two  top  valence  bands  which,  as  a  result, 
are  parabolic.  We  adopt  an  isotropic  phonon  spectrum  with  spher¬ 
ically  averaged  phonon  velocities  for  the  three  polarization 
branches.  Tiersten  and  Lawaetz  do  use  full  solutions  of  the 
dynamical  matrix  for  the  phonon  spectra  and  polarization  vectors. 
Lawaetz  shows,  however,  the  results  he  obtained  are  insensitive 
to  this  detail*4*. 

In  Figs.  1  and  2  we  show  the  solution  for  that  part  of 
the  distribu*  an  function  affected  by  the  application  of  electric 
field  only,  5,  for  the  L  and  H  bands  of  Ge .  It  is  convenient 
to  display  9,  Eq .  (29),  where 

0NX(E)aE  [4,MX-<E,r1.GMX-IE>  •  l7'1> 

MX' 

and  has  the  units  of  (Ry  a.u.).  In  practical  units  Eq.  (30) 
is 
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Fig.  2.  Cubic  harmonic  decomposition  of  the  perturbed  part  0 
of  the  distribution  function,  Eq.  (74) ,  as  a  function 
of  hole  energy,  for  the  heavy  hole  band  of  germanium. 
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where  t  is  the  numerical  value  of  the  absolute  temperature 

T  =  (t)  ( °  K )  ,  (76a) 

and  9  is  the  numerical  value  of  9, 

9  =  6  (Ry  a. u. )  .  (76b) 

Since  f^(k)  is  dimensionless,  5  has  the  units  of  inverse  electric 
field. 

From  Figs.  1  and  2  we  see  that  t  is  fairly  well  converged 
by  the  time  the  last,  5^,  cubic  harmonic  is  taken  into  account. 

9  for  the  heavy  hole  band  is  remarkably  constant  over  the  entire 
energy  rj.nge,  whereas  0  for  the  light  hole  band  has  an  appreciable 
energy  dependence  at  very  low  energies.  This  appears  to  be  due 
to  the  closer  proximity  of  the  L  band  to  the  S  band  resulting  in 
a  greater  nonparabol ici ty  for  the  L  band  than  for  the  H  band.  In 
practice,  the  H  band,  being  more  flat  than  the  L  band,  contains 
more  holes  at  any  given  temperature  so  that  the  transport  proper¬ 
ties  of  Ge  (and  Si)  are  greatly  influenced  by  the  less  mobile 
holes  in  the  H  band.  Yet,  Fig.  1  does  indicate  that  even  in  Ge 
one  cannot  treat  the  L  band  as  being  parabolic.  In  fact,  one  can 
consult  Fig.  1  of  Kane's  paper ^ 40 '  to  ascertain  that  the  H  band  is 
accurately  parabolic,  while  the  L  band  is  not.  From  the  dimen¬ 
sional  analysis,  Eq.  (70)  ,  we  expect  for  parabolic  bands  that 
^constant.  In  fact,  this  is  closely  reflected  in  Fig.  2  for  the 
H  band,  except  at  the  very  lowest  of  energies  where  there  may 

be  difficulties  due  to  the  neglect  of  the  phonon  energy  in  the 

-h 

argument  of  the  delta  function,  Eq.  (3).  Since  (72),  any 

strong  energy  dependence  in  4> ,  or  S  and  H,  will  reflect  itself  in 
the  x  vs*  E  curves,  which  we  show  next. 

In  Figs.  3  and  4  we  show  the  solution  for  the  part  of  the 
distribution  function  affected  by  the  application  of  the  magnetic 
field,  x>  f°r  the  L  and  H  bands  of  Ge .  From  Eq.  (43)  we  define 
an  auxiliary  quantity,  £,  through 
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E(0.0045eV) 


Fig.  3.  Cubic  harmonic  decomposition  of  the  perturbed  part  £ 
of  the  distribution  function,  Eq.  (77) ,  as  a  function 
of  hole  energy,  for  the  light  hole  ’-and  of  germanium. 
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Fig.  4.  Cubic  harmonic  decomposition  of  the  perturbed  part  £ 
of  the  distribution  function,  Eq .  (77),  as  a  function 
of  hole  energy,  for  the  heavy  hole  band  of  germanium. 


5v,r.  .  /  4 tt  \ j  4Tre2  r<5  ~  ~  . 

XNX(E)  “  (kBT  j  fi  ENX(E)  (FxB)  *ey 


(77) 


2  3 

whe-e  E„\(f)  is  in  (Ry  a.u.  ).  In  practical  units 


o v  r  1.6 2 4x10 9  cm3  ^6  ,  r  "  ~  ~  i 

<nx(E>  -  — "7 - 3 -  '-nx(E)  [(FXB)  -e  l 

t  V  sec 


(78) 


where  t  is  defined  in  (76a),  and  (Ry2a.u.3).  Since  the 

Vsec  V  —l 

B  field  is  in  - T  and  F  ^  v  has  the  units  of  (FB)  to 

cm2  cm 

make  f„  dimensionless 
N 

Figs.  3  and  4  mostly  meet  our  expectations  of  displaying 

_L 

the  E  behavior,  Eq.  (72).  The  part  of  E  for  the  L  band 
displays  a  rather  violent  downturn  near  the  origin.  This 
phenomenon  can  be  traced  to  inflections  in  the  corresponding 
6  curve.  Fig.  1.  we  would  like  to  point  out  that  this  behavior 
occurs  over  a  very  limited  range  of  energies  near  the  origin 
and  should  not  affect  our  results  for  but  the  lowest  temperatures, 
T<10K  or  so.  In  part,  we  can  attribute  the  behavior  to  the 
known  shortcomings  of  the  theory  for  the  lowest  of  energies  where 
some  assumptions  of  the  theory  break  down.  of  these,  the  most 
important  one  is  the  assumption  of  the  elastic  scattering 
between  holes  and  phonons,  Eq.  (3),  for  low  hole  energies. 

To  further  assure  ourselves  of  the  validity  of  our  calcu¬ 
lation,  and  to  gain  a  further  insight  into  its  results,  we  present 
another  intermediate  result  before  displaying  mobility  curves  for 
Ge.  The  integrands,  apart  from  the  Boltzmann  factor,  in  Eqs . 

(53a  and  b)  are  known  functions  of  energy  in  the  parabolic  limit. 
We  define 


VE> 

.  E 

X 

9nx<e>5nU<e» 

in  (Ry  a.u.-1) 

(79) 

dn(e) 

=  4fr<a 

o 

1  rN 1 

in  (Ry-1a.u.-3)  , 

(80) 
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and 


Wn(E)  =  Sc^(E)g^(E)  in  (Ry2  a.u.).  (81) 

D^CE)  is  related  to  the  density  of  states  and  for  parabolic  bands 

D  ^  r  %  yh  5F  %  '  (82) 

whereas 

U  'v.  6g  v  (const.)  x  E  v  E  ,  (83) 

W  ^  Eg  ^  E^E  v  E1'2  .  (84) 

Figs.  5,  6,  and  7  display  U,  D,  and  W,  respectively,  for  the 
light  and  heavy  hole  bands  of  germanium.  In  spite  of  the  com¬ 
plicated  0  dependence,  displayed  in  Fig.  1,  the  U  curve  for  the 
L  band  is  close  to  being  linear  in  E  except  for  a  small  inflection 
for  E<0.0045  eV.  The  plot  for  the  H  band  is  almost  indistinguish¬ 
able  from  a  straight  line.  As  expected  U  >UT ,  by  about  a  factor 
of  four,  since  the  heavy  hole  band  has  a  greater  density  of  states 
than  the  L  band.  Since  ,  the  conductivity  mobility  will 

be  influenced  by  the  H  band  more  than  by  the  L  band.  It  turns  out 
that  the  effect  of  the  \  sum  and  the  multiplication  by  g  in  Eq. 

(79)  is  to  smooth  out  the  small  imperfections  in  the  calculation 
of  Figs.  1  and  2,  at  small  energies. 

Fig.  6  displays  the  density-of-states  like  curves,  , 
for  the  L  and  H  bands  of  Ge.  As  expected  these  do  behave  as  /E. 
The  densities  of  states  are  in  the  approximate  ratio  DtJ/DT'v20,  so 
that  there  are  many  more  H  holes  than  L  holes  at  any  given  temper¬ 
ature.  Fig.  7  for  W^,  which  enters  the  Hall  mobility  calculation, 
again  displays  the  expected  /F  behavior  and  WT/Wtl'v2,  therefore, 
the  Hall  mobility  will  be  dominated  by  the  L  band.  In  all  the 
cases,  Figs.  5-7,  the  small  energy  dependence  was  verified  and 
holds  true  to  good  accuracy.  This  is  in  spite  of  rather  compli¬ 
cated  behavior  of  9  and  E  in  Figs.  1-4.  In  part,  this  is  due  to 
the  averaging  effect  in  summing  over  \  in  Eqs .  (79)  and  (81), 

and  to  the  strong  band  effects  contained  in  g(E).  This  does  not 
mean,  of  course,  that  any  part  of  the  treatment  can  or  should  be 
simplified  since  every  step  leading  up  to  the  final  result  has  a 
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Fig.  5. 


Temperature  independent  part  of  the  integrand  in  the 
numerator  of  Eq.  (53a)  for  u  ,  defined  by  Eg.  (79). 
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f 


Fig.  6.  Temperature  independent  part  of  the  integrand  in  the 
denominators  of  Eqs .  (53a  and  b)  for  u  and  y„ , 
defined  by  Eg.  (80) . 
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Fig.  7.  Temperature  independent  part  of  the  integrand  in  the 
numerator  of  Eg.  (53b)  for  y  ,  defined  by  Eq.  (81). 

H 
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quantitative  bearing  on  its  accuracy.  Finally,  Figs.  5-7  validate 
our  low  energy  interpolation  procedure  discussed  in  Eqs.  (69)— (73) . 


Figs.  8  and  9  show  the  total  and  partial  conductivity  and 
Hall  mobilities,  respectively,  as  a  function  of  absolute  temper¬ 
ature  for  Ge .  The  partial  mobilities  are  defined  by  the  following 
equation 


U 


c,H 


L 

_N _ 

E 


PNK 


N 

c ,  H 


P 


N 


(85) 


N 

where  u  „  are  the  partial  mobilities.  Eq.  (85)  is  motivated  by 

C  /  rl 

a  model  where  the  bands  conduct  current  as  a  parallel  resistor 
network.  In  this  model,  therefore,  the  band  conductivities  are 
additive  as  though  the  bands  were  noninteracting.  Since  we  have 
shown  that  the  interband  scattering  is  as  important  as  the  intra¬ 
band  scattering,  the  motivation  for  the  model  does  not  hold  any 
longer.  Yet,  we  can  still  define  the  partial  band  mobilities 
through  Eq .  (85)  and  draw  some  conclusions  from  the  results.  In 

the  case  of  both  mobilities  the  light  hole  band  mobility  is  the 
larger  of  the  two.  Figs.  8  and  9,  with  u  /u  ^5  and  u„/|J„^40. 

C  C  H  H 

Physically,  the  larger  light  hole  band  mobilities  are  the  conse- 
quences  of  its  having  a  smaller  effective  mass  ,  Mu/M.,.^8.  In 
the  case  of  the  Hall  mobility  the  L  band  becomes  the  dominant  band 
in  spite  of  the  fact  that  it  contains  twenty  times  fewer  holes 
than  the  H  band. 


Let  us  now  compare  the  results  of  our  theoretical  formula¬ 
tion  of  the  problem  and  calculational  formalism  with  the  results 
of  calculations  by  Tiersten  and  Lawaetz.  Tiersten  calculated  the 
conductivity  mobility  for  Ge  in  the  parabolic  band  approximation 
using  the  H  and  L  bands  only.  As  we  have  shown,  the  L  band's  non- 
parabolicity  has  its  strongest  effect  on  the  Hall  mobility.  Using 
the  parameters  listed  in  Table  1  Tiersten  obtained  the  expected 
T  1-5  temperature  dependence  for  the  conductivity  mobility,  as 
did  Lawaetz,  with 

UC(100K)  =  3.13  x  10^  cm^/Vsec. 
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Fig.  8.  Calculated  total  and  partial  conductivity  mobilities, 
in  the  acoustic  phonon  limited  regime,  as  a  function 
of  temperature. 
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Fig.  9.  Calculated  total  and  partial  Hall  mobilities,  in  the 
acoustic  phonon  limited  regime,  as  a  function  of 
temperature . 


Lawaetz  using  essentially  the  same  formalism  as  ours,  but  with 
two  parabolic  bands,  states  that  his  100K  value  for  pc  differs 
from  Tiersten's  value  by  about  2%.  Our  value  for  the  acoustic 
phonon  iimited  conductivity  mobility  at  this  temperature  is 

yc(100K)  =  2.89  x  10^  cm^/Vsec  , 

which  differs  from  Tiersten's  result  by  about  7.7%. 

This  small  difference  can  be  accounted  for  by  several 
factors  which  distinguish  our  respective  calculations.  One  of 
these  might  be  the  precision  with  which  the  transition  proba¬ 
bilities  were  fitted  to  the  double  cubic  harmonic  series.  We 
use  a  far  larger  sampling  of  (k,k')  scattering  directions  than 
Tiersten  and  a  larger  number  of  cubic  harmonic  pairs  for  the 
fit  .  In  addition,  we  do  not  treat  the  L  band  as  parabolic 
which  is  a  progressively  worse  approximation  as  higher  galvano- 
magnetic  coefficients  are  calculated.  On  the  other  hand, 
Tiersten's  more  rigorous  treatment  of  the  phonon  spectrum  is  an 
improvement  over  our  isotropic  phonon  spectrum  approximation. 
However,  we  believe  that  it  is  the  parabolic  L  band  approxima¬ 
tion  which  accounts  for  the  bulk  of  the  difference. 

Fig.  10  displays  both  the  calculated  r-factor  for  germanium, 
data  to  be  discussed  later,  and  exponents  a  , a„  in 

C  ri 

-a 

uc  =  AT  C  (86) 

— r i 

UH  =  BT  H  .  (87) 

The  temperature  exponent  of  the  conductivity  mobility,  a  ,  starts 
out  at  10K  with  ac  =  1.498  and  rises  monotonically  to  ac  =  1.575 
at  120K.  The  temperature  exponent  for  the  Hall  mobility, 
rises  steadily  from  a„  =  1.420  at  10K  to  ct„  =  1.545  at  120K. 

ri  ri 

This  represents  a  rather  small  temperature  variation  about  the 
value  a  „  =  1.5  for  strictly  parabolic  bands. 

C  t  ri 

It  is  more  difficult  to  compare  our  results  for  uH  and  r 
with  those  of  Lawaetz  (Tiersten  did  not  calculate  these  quanti¬ 
ties)  .  Lawaetz  does  obtain  a„  =  1.5  which  is  the  correct 

H 
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Fig.  10.  Top  scale:  Calculated  temperature  exponents  a 
for  the  conductivity  mobility,  and  aH  for  the 
Hall  mobility  (solid  line) . 

Lower  scale:  The  calculated  r-factor  (solid  line) 
and  experimental  data  of  Beer  and  Willardson- 
triangles  (Ref  53,  impurity  content  IxlO3 3cm-3 ) , 
Goldberg,  et  al. -circles  (Ref.  54,  sample  D, 
impurity  content  2 . 4x10 1 3cm“3 ) . 
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exponent  for  both  mobilities  in  the  parabolic  band  regime.  Our 
aH  is  lower  than  1.5  for  T<70K  thus  reflecting  the  dominant 
influence  of  the  light  hole  band  to  u„.  Lawaetz  provides  the 
value  of  r  =  1.60  for  the  choice  of  the  deformation  potential 
parameters  b  =  -2eV,  d  =  -7eV,  and  a  value  of  the  dilatational 
deformation  potential  a>0  which  will  yield  the  lowest  value  for 
r  for  acoustic  mode  scattering  only.  The  value  of  r  =  1.60, 
therefore,  represents  the  lowest  value  Lawaetz  calculated  for 
acoustic  phonon  limited  mobility  with  fixed  b  and  d,  and  a 
positive  value  for  a.  From  Fig.  3  of  Ref.  2  we  estimate  that 
for  the  parameters  in  Table  1  Lawaetz  obtained  r^l.6-1.7. 

The  experimental  value  for  r  that  Lawaetz  used  for  compar¬ 
ison  with  the  theory  is  the  value  of  r  =  1.36  +0.07  at  81K 
obtained  by  Beer  and  Willardson ^ 5 3 '  (sample  impurity  content 
lxl013cm  3).  This  is  far  below  Lawaetz's  r^l.6-1.7.  On  the 
other  hand,  our  result  at  81K  is  r  =  1.43,  Fig.  10,  which  is 
very  close  to  the  r  of  Beer  and  Willardson,  with  the  parameters 
used  by  Lawaetz. 

We  can  suggest  primarily  two  reasons  for  the  different 
results  of  the  two  theoretical  calculations.  As  we  have  shown, 
the  light  hole  band  is  quite  nonparabolic  and  dominates  the  Hall 
mobility.  The  L  band,  therefore,  alters  uu  and  r  from  the  values 
calculated  in  the  parabolic  approximation.  We  also  have  one  mis¬ 
giving  about  Lawaetz 1  a  chosen  form  for  his  scattering  matrix  K 
(S  in  our  notation).  His  K  matrix,  Eq.  3.15  of  Ref.  4,  uses  the 
completeness  relation  for  cubic  harmonics  in  the  second  term  of 
K  involving  the  sum  over  uv.  Formally,  this  is  correct,  but 
computationally  one  has  to  achieve  convergence  in  the  sum  over  an 
additional  angular  momentum  index.  The  sum  over  u  can  be  elimi¬ 
nated  by  using  the  completeness  relation  so  there  will  be  only 
one  sum  over  v  to  carry  to  convergence.  We  feel  this  would  be  a 

safer  procedure  considering  that  in  order  to  calculate  one 

- 1  .  ” 

must  use  S  twice,  Eqs .  (29)  and  (43). 
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Fig.  10  also  displays  the  data  for  sample  D  of  Goldberg, 
et  al.*54^  with  the  impurity  concentration  of  2.4xl013cm~3  (their 
sample  G  was  less  pure  than  sample  D) .  Older  data  of  Morin(55) 
for  r  is  not  shown,  since  the  impurity  content  is  unknown,  but  it 
generally  falls  below  the  data  of  Goldberg,  et  al.  Thus,  by 
comparison  with  the  data  of  Goldberg,  et  al.  the  calculated  r 
factor  is  about  12%  too  high.  It  appears  that  the  data  of  Beer 
and  Willardson  coming  from  a  purer  sample  may  be  more  representa¬ 
tive  of  a  high  purity  sample  with  less  of  ionized  impurity  scatter¬ 
ing.  By  contrast,  Nagakawa  and  Zukotynski (  10 '  calculate  r^2 . 0  in 
this  temperature  range,  far  above  the  measured  values.  Ref.  10 
also  shows  that  ionized  impurity  scattering  has  the  effect  of 
dramatically  lowering  the  value  of  r. 

Although  it  has  not  been  our  aim  to  provide  a  critical  com¬ 
parison  between  experiment  and  theory  for  Ge ,  we  can  offer  several 
suggestions  concerning  the  results.  First  of  all,  the  data  used 
for  comparison  in  Fig.  10  is  in  the  temperature  range  where  the 
optical  phonon  scattering  should  come  in  for  germanium.  There¬ 
fore,  the  data  may  not  be  representative  of  the  temperature  regime 
for  the  acoustic  phonon  limited  mobility  T<80K.  Second,  we  would 
have  liked  to  use  for  comparison  data  obtained  from  higher  purity 
samples,  on  the  order  of  1011cm-3,  for  low  temperatures,  T<80K. 

No  such  data  for  the  same  sample  seems  to  exist  for  Ge  where  both 
the  high  and  low  magnetic  field  limits  were  explored  to  obtain 
the  conductivity  and  Hall  mobilities,  respectively.  Third,  the 
values  of  the  deformation  potential  parameters  for  Ge  are  still 

in  doubt.  In  particular,  the  compilation  of  these  parameters  by 

(28) 

Wiley  suggests  that  our  value  for  d  is  probably  too  high. 

We,  however,  do  not  feel  that  the  experimental  situation  here  has 
been  resolved  enough  to  warrant  another  calculation  for  Ge  at 
this  time. 

Fig.  11  presents  the  comparison  between  our  calculated 
conductivity  mobility  and  results  of  measurements  by  Ottaviani, 
et  al.(56)  and  Brown  and  Bray<57^.  A  rather  close  agreement  is 
found  for  T<80°K,  while  for  T>80K  one  can  see  already  the 
influence  of  the  optical  phonon  scattering  in  the  data  of  Brown 

(54) 

and  Bray.  Also,  the  20°K  uc  point  of  Ottaviani  et  al.  may 
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Fig.  11.  Calculated  conductivity  and  Hall  mobilities  (solid 
lines)  and  the  data  for  uc  of  Ottaviani,  et  al.- 
c^rc^es  (Ref.  56,  time-of -f 1 iqht  method,  impurity 
content  <lxl011cm~3)  r  and  Brown  and  Bray-triangles 
(Ref.  57,  impurity  content  1 . 25x10 1 3cm"^ ) . 
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be  indicating  the  presence  of  ionized  impurity  scattering. 

Owing  to  the  lack  of  the  y^  and  uc  data  on  the  same  high  purity 
sample  for  low  temperatures,  we  cannot  present  a  comparison  with 
theory  for  y^.  For  this  purpose  it  is  best  to  consult  Fig.  10 
for  the  r-factor. 

The  calculation  for  germanium  has  achieved  its  purpose  as 
a  testing  ground  for  the  theoretical  and  computational  procedures 
We  have  reproduced  the  results  of  earlier  yc  calculations  of 
Tiersten  and  Lawaetz  and  improved  considerably  the  results  for 
yH  by  including  the  nonparabolicitv  for  the  L  band.  The  degree 
of  quantitatitive  agreement  this  calculation  has  achieved  with 
experiments  does  indicate  that  a  first-principles  approach  can 
be  quantitatively  successful  and  confirms  the  fundamental  correct 
ness  of  the  deformation  potential  theory.  This  viewpoint  will  be 
reinforced  in  the  next  section  on  Si  where  a  more  reliable  set  of 
deformation  potential  parameters  will  be  used  and  comparison  made 
with  our  experimental  results. 
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SECTION  V 

RESULTS  FOR  SILICON 


1.  EXPERIMENT 

The  material  used  for  this  experiment  was  all  ultra-high 

purity,  vacuum  float  zoned,  p-type  silicon.  The  residual  acceptor 

in  this  material  is  boron,  which  was  found  by  conventional  Hall 

12  —  3 

analysis  to  be  present  in  concentrations  in  the  range  of  10  cm 

or  less.  The  total  residual  donor  concentration  was  in  the  low 

1011cm  3  range  for  all  but  one  of  the  samples  studied.  The  donors 

are  most  likely  phosphorus  but  were  not  identified.  Both  Hall  bars 

and  van  der  Pauw  cloverleaf  samples  were  cut  from  wafers  of  several 

boules.  The  Hall  coefficients  of  these  samples  were  measured  over 

the  magnetic  field  range  of  0.003  to  2.0  T  at  several  temperatures 

between  20  and  50K.  The  typical  range  of  ojct  at  20K,  determined 

by  assuming  cocx  =  vcH,  was  between  0.05  and  30,  so  we  were  clearly 

able  to  reach  both  the  high  and  the  low  field  limits  with  our 

magnet.  Accurate  measurements  of  Ru(°°)  were  prevented  for  temper- 

atures  above  50K  by  the  decrease  in  the  mobility,  which  forced  the 

high  field  values  of  wct  considerably  below  10.  Because  of  some 

lack  of  certainty,  even  at  the  lowest  temperatures,  as  to  the 

precise  value  of  RH(°°),  the  high  field  data  at  each  temperature 

were  fitted  to  a  simple  form  to  yield  the  saturation  value.  It  is 

well  known  that  the  Hall  coefficient  must  be  an  even  function  of  H, 

_  2 

which  means  that  R  (°°)  must  be  approached  as  H  .  We  therefore 

11 

fitted  our  high  field  data  to  the  form 

Rh(H)  =  Rh(°°)  +  AH'2  (88) 

with  rh(°°)  and  A  as  adjustable  parameters.  Only  values  of  R^ 
for  coct>4  were  used  at  any  temperature  for  this  fitting.  The  low 
temperature  limit  for  the  measurements  was  determined  by  the  input 
impedance  of  our  Hall  effect  system  which  is  about  1014ft.  This 
limits  the  sample  resistance  to  10  ft  and  the  sample  resistivity 
to  about  101Oft-cm.  This  system  is  described  elsewhere ( 58 * . 


Because  most  Hall  effect  measurements  are  made  in  the  low 

/  c 

field  limit,  to  t<<1,  where  R„  also  approaches  a  constant  value 

C  ri 


48 


we  will  define  a  magnetic  field  independent  r  for  the  purposes 
of  this  report  as  the  ratio  of  Ft,  in  the  two  limits  or 


y°> 

V"> 


Rh  (  0)  ne 


(89) 


where  RH(0)  and  RH(°°)  are  respectively  the  low  and  high  field 
limits  of  the  Hall  coefficient.  In  a  similar  manner,  the  Hall 
mobility  will  be  considered  as  the  low  field  Hall  mobility. 

Fig.  12  shows  one  of  our  cleaner  field  sweeps  plotted  as 
RH  (“ct)  /RH  (oo)  versus  ojct  .  In  this  form  the  low  field  plateau 
represents  r.  The  solid  line  is  drawn  to  show  the  trend  in  the 
data  and  does  not  represent  any  attempt  to  fit  the  data  to  a 
functional  form  throughout  the  entire  field  range.  Fig.  13,  which 
is  a  collection  of  all  such  lines  for  the  sample  shown  in  Fig.  12, 
clearly  shows  a  temperature  dependence  in  r.  This  temperature 
variation  was  seen  in  all  the  samples  we  studied.  Fig.  14  shows 
the  temperature  variation  of  r  for  our  four  samples.  The  open 
symbols  are  Hall  bars  and  the  solid  ones  are  van  der  Pauw  clover- 
leaves.  All  four  samples  came  from  different  boules  but  were 
all  of  similar  purity.  Table  2  gives  the  results  of  conventional 
(r=l)  Hall  analysis  of  these  samples. 


TABLE  2 

HALL  ANALYSES  OF  SAMPLES 


Sample 

Shape 

NB(cm'3) 

Np (cm  3) 

Vnd 

1300-V 

v.d.P. 

9.14X1011 

3 . 33X1011 

2.74 

1202-H 

Hall 

6. 57X1011 

3. 96xlOU 

1.66 

2001-V 

v.d.P. 

2.62xl012 

2.08xl012 

1.26 

4802-H 

Hall 

5.93X1011 

2.36X1011 

2.54 

Although  the  scatter  in  Fig.  14  is  considerable,  the 
difference  between  the  Hall  bars  and  the  van  der  Pauw  samples  could 
be  real.  The  Hall  bar  is  a  more  satisfactory  configuration  for 
high  magnetic  field  measurements  than  the  cloverleaf  shape  because 
the  current  path  is  less  tortuous  and  therefore  less  likely  to  be 
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distorted  in  large  fields.  Some  aspects  of  probe  effects  and 
sample  shape  are  discussed  by  Frederikse  and  Hosier*60^. 

The  dominant  feature  in  the  curves  in  Fig.  13  is  the  hump 
near  gjct=1.  The  position  of  the  hump  shows  little  variation  with 
temperature  or  between  samples  but  the  height  does  tend  to  vary 
between  samples,  those  with  the  larger  r  at  a  given  temperature 
show  a  slightly  higher  hump.  This  feature  somewhat  resembles  the 
structure  in  the  curves  of  r  versus  H  for  germanium  calculated  by 
Beer  and  Willardson ( 6 1 }  which  is  a  result  of  the  multiple  nature 
of  the  valence  bands.  Since  silicon  and  germanium  have  similar 
multiple  valence  bands,  we  feel  that  this  is  the  likely  source  of 
this  feature.  Unfortunately,  the  scatter  in  our  data  is  too  large 
to  make  any  quantitative  statements  about  this  interesting  structure 

During  the  conventional  Hall  effect  versus  temperature 
measurements  made  on  these  samples  for  analysis,  we  also  measured 
the  Hall  mobility  from  20  to  380K  at  0.1T  which  is  shown  in 
Fig.  15  for  one  sample.  At  the  low  temperature  end  the  mobility 

_CtH 

appears  to  obey  a  simple  power  law,  u  “T  ,  where  a„  is  a  positive 

n  n 

number.  This  dependence  holds  up  to  about  100K  where  there  is  a 
transition  to  a  steeper  negative  slope  in  the  log  uH  versus  log  T 
plots.  This  break  is  usually  attributed  to  the  onset  of  optical 
phonon  scattering  while  the  shallower  slope  is  characteristic  of 
acoustic  phonon  scattering.  There  is  no  indication  of  a  turnover 
toward  a  positive  slope  at  the  low  end  which  would  indicate  ionize  I 
impurity  scattering.  This  is  consistent  with  the  fact  that  there 

12  -  3 

are  only  10  cm  ionized  centers  in  this  material.  The  concentra¬ 
tions  of  neutral  scattering  centers  are  in  the  1015-1016cm-3  range 
and  are  primarily  non-active  oxygen  and  carbon  impurities  that  were 
measured  in  similar  material  by  optical  absorption.  The  scattering 
in  the  temperature  range  of  this  experiment  is  therefore  exclusively 
due  to  acoustic  phonons. 
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2. 


COMPARISON  WITH  PREVIOUS  THEORETICAL  RESULTS 


Many  attempts  have  been  made  over  the  years  to  calculate  the 
Hall  factor  and  the  conductivity  mobility  for  acoustic  phonon 
scattering.  The  simplest  calculation  is  a  Drude-type  approach 
which  assumes  a  constant,  isotropic  scattering  time  and  free 
electrons,  that  is,  a  single,  spherical,  parabolic  band.  This 
approach  yields  tJ a  simple  result  that  r=l  for  all  magnetic  fields. 
This,  of  course,  implies  that  the  conductivity  and  Hall  mobilities 
are  identical.  The  temperature  dependence  of  the  mobility  in  this 
model  is  due  to  the  statistics  of  the  phonons  and  electrical 
carriers.  It  is  only  when  the  scattering  time  is  a  function  of 
energy  or  velocity  that  r  deviates  from  unity.  The  standard  text¬ 
book  value  of  3tt/8  for  acoustic  phonon  scattering  is  obtained  by 

•  ,  —  \c 

assuming  the  scattering  time,  t,  is  proportional  to  E  ,  again 
with  a  spherical  and  a  parabolic  band^62*.  Lax  and  Mavroides ( 6 3 ' 64 5 
included  the  warpage  of  the  valence  band;  i.e.,  nonspherical 
energy  surfaces,  but  did  not  include  nonparabolicity ,  which  is 
important  for  the  valence  bands  in  this  material.  They  used  a 
relaxation  time  of  the  form  t=£E  \  The  value  of  r  they  obtained 
was  again  independent  of  temperature.  More  recently,  Nakagawa  and 
Zukotynski * 10) ,  using  a  relaxation  time  approximation,  calculated 
the  conductivity  mobility  of  p-type  silicon  for  acoustic  phonon 
scattering.  They  included  both  the  warpage  and  the  nonparabolicity 
of  the  two  major  valence  bands.  In  addition,  they  used  a  more 
realistic  scattering  time. 

We  used  (°°)  and  the  zero  field  conductivity  to  determine 
the  conductivity  mobility,  and  in  a  similar  manner  determined  the 
true  low  field  Hall  mobility  using  R^CO).  Fig.  16  shows  these  two 
mobilities  for  the  sample  presented  in  Figs.  12  and  13.  The 
temperature  dependences  are  also  given  assuming  a  simple  exponential. 
The  exponent  of  T  in  the  Hall  mobility  for  all  samples  studied 
varied  between  1.6  and  1.7  while  that  for  the  conductivity 
mobilities  varied  between  1.7  and  1.8.  For  a  given  sample  the 
exponent  for  the  Hall  mobility  was  always  the  smaller  positive 
number  of  the  two. 

i 

i 
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Fig.  16.  Mobilities  of  sample  1202-H.  Open  circles: 

conductivity  mobilities;  solid  circles:  Hall 
mobilities . 
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From  the  simple  physical  arguments  the  temperature  depen¬ 
dence  of  the  Hall  factor  can  be  due  to  at  least  one  of  two  possible 
effects.  Firstly,  the  nonparabolicity  of  the  valence  bands  com¬ 
bined  with  the  warpage  will  produce  a  temperature  variation  and 
this  can  be  seen  quite  simply.  As  shown  by  Lax  and  Mavroides, 
warpage  alone  produces  a  deviation  from  one  in  the  Hall  factor  but 
if  the  bands  are  parabolic,  there  will  be  no  change  with  temper¬ 
ature  because  the  shape  remains  constant  as  the  energy  increases. 
If,  however,  there  is  a  deviation  from  simple  parabolicity ,  the 
shape  of  the  constant  energy  surfaces  will  vary  with  increasing 
energy.  As  the  temperature  increases,  more  and  more  carriers  will 
populate  progressively  deeper  energy  states  and  as  they  do  they 
will  see  differently  shaped  energy  surfaces,  thus  changing  the 
overall  Hall  factor.  The  other  possible  source  of  the  temperature 
variation  could  be  a  change  in  the  fractional  density  of  carriers 
in  the  three  valence  bands.  The  light  and  heavy  hole  and  the 
split-off  bands  each  have  different  warpages  and  nonparabolicities 
so  each  taken  individually  will  have  a  different  r.  If  the  dis¬ 
tribution  of  carriers  amongst  these  bands  varies  with  temperature, 
the  total  r  will  therefore  also  have  to  vary. 

Assuming  the  simple  form  r=aT  ,  we  obtained  a  fit  to  the 

measured  r's  for  the  two  Hall  bar  samples.  The  result  of  this  fit 

is  the  solid  line  in  Fig.  14.  Our  best  fit  gave  a=0.85  and 

b=0.082.  We  used  this  simple  fit  to  correct  an  anomaly  that  has 

been  reported  in  the  p  vs.  1000/T  data  for  ultra  pure,  p-type 

silicon  and  tentatively  attributed  to  the  effect  of  neglecting 

the  Hall  factor  in  the  data  analysis.  Baron,  et  al.*65)  reported 

that  the  saturation  region  in  this  material  does  not  behave  as 

expected.  Instead  of  remaining  flat  after  saturation  is  achieved, 

the  hole  concentration  sometimes  decreases  slightly  just  before  the 

steep  rise  due  to  band  gap  excitation  and  residual  deep  centers 

sets  in.  This  occurs  roughly  between  60  and  100K  for  samples  with 

12-3 

boron  concentrations  below  10  cm  .  This  temperature  range  is 
outside  that  covered  by  our  measurements  but  is  still  mostly  in 
the  region  dominated  by  acoustic  phonon  scattering  so  it  is  logical 
to  assume  that  the  functional  dependence  of  r(T)  should  not  change 
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drastically  from  that  which  we  measured  at  the  lower  temperatures. 
Following  this  assumption,  we  took  the  p  vs.  1000/T  data  from 
the  conventional  Hall  effect  run  on  sample  1202-H  and  corrected 
the  p's  with  our  empirical  r  at  each  temperature.  The  results  are 
shown  in  Fig.  17,  which  is  a  greatly  expanded  plot  of  the  region 
of  interest.  The  inclusion  of  the  Hall  factor  has  removed  the 
anomaly.  In  addition  to  this,  we  corrected  all  the  data  for  this 
run  and  then  refitted  the  data  to  determine  the  boron  concentration 
(the  effective  mass  of  Ref.  48  was  used  with  the  charge  balance 
equation  to  fit  this  data) .  The  inclusion  of  the  Hall  factor 
changed  the  boron  concentration  and  energy  from  6.57  x  1011cm  3 
and  .0448eV  to  7.09  x  1011cm-3  and  .0449  respectively.  The  cal¬ 
culated  donor  concentration  remained  unchanged  at  3.96  x  10i:Lcm  3. 
The  largest  change  is  in  the  boron  concentration  which  has 
increased  by  about  10  percent. 


3.  PRESENT  THEORETICAL  RESULTS  FOR  SILICON 

The  mobility  calculation  for  silicon  has  been  performed 

using  the  parameters  listed  in  Table  3.  The  deformation  potential 

parameters  given  in  Table  3  are  not  the  ones  used  in  the  two  pre- 

(4345) 

vious  publications  '  .  The  two  sets  of  numbers  come  from  two 

(28  52) 

different  compilations  by  Wiley  '  .  The  set  used  here  is  the 

(28) 

more  recent  one  and,  moreover,  is  in  good  agreement  with  the 

latest  piezospectroscopic  measurements  on  Si130'38*. 


TABLE  3 

EXPERIMENTAL  INPUT  PARAMETERS  USED  IN 
THE  MOBILITY  CALCULATION  FOR  SILICON 


Valence  Band 
Parameters3 
(dimensionless) 


Deformation 
Potentials*5  (eV) 


Phonon  Parametegs 
dc2  (dyne/cm2 ) 


A 

-4.27 

a 

2.1 

B 

-0.63 

b 

-2.2 

C 

4.93 

d 

-5.3 

s=transverse  e.BOxlO"*"1 
s=longitudinal  18.852x10 


3  Ref.  10 

b  Ref.  28,  p=mass  density,  cs=phonon  rueed  for  branch  s 
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1000/T  (K'1) 

Hole  concentration  vs.  1000/T  for  sample  1202-H. 
Solid  circles:  uncorrected  data  (left  axis) ;  open 
circles:  same  data  corrected  for  Hall  factor 
(right  axis) . 


Fig.  18  provides  functions  0,  Eq.  (74)  for  the  L  and  H 
bands  of  silicon.  As  was  the  case  in  germanium,  the  6^  components 
of  9  are  the  largest  and  the  6^  and  <S£  components  the  smallest. 
Since  the  signs  of  the  6^  contributions  are  positive  the  distri¬ 
bution  functions  for  both  bands  are  displaced  in  k  space  in  the 
direction  of  the  F  field.  6's  for  the  H  band  are  somewhat  smoother 
than  the  0's  for  the  L  band  indicating  that  the  H  band  is  the  more 
parabolic  of  the  two.  The  departure  of  the  0  curves  from  a  B^con- 
stant  behavior  can  be  attributed,  in  large  part,  to  the  fact  that 
the  expansion  coefficients  A,  Eq.  (3)  are  energy  dependent  them¬ 
selves.  Clearly,  the  closer  proximity  of  the  L  band  to  the  5  band 
is  the  cause  of  the  greater  nonparabolicity  of  the  L  band. 

Figs.  19  and  20  display  the  energy  dependence  of  the  E,  parts 
of  the  distribution  functions  for  the  L  and  H  bands,  respectively. 
Overall,  the  E  contributions  do  behave  in  the  approximate  E  2 
manner.  The  6'  part  of  E  for  the  H  band  starts  out  at  first  like 
-(E)  and  then  reverts  to  the  (E)  ^  behavior.  The  same  trend  was 
seen  in  the  6^  part  of  E  for  the  L  band  of  Ge .  We  see  that  the  E 
functions  are  not  as  well  converged  as  in  germanium.  As  far  as 
mobility  calculations  are  concerned  the  E  convergence  in  X  is  not 
the  important  question.  What  must  be  converged  is  the  sum  of 
products 

Z4x<E,6m(E)  -  <90> 

X 

for  the  conductivity  mobility,  Eq.  (53a) ,  and 

£%“(E)4xIE>  ,  (9i) 

for  the  Hall  mobility,  Eq.  (53b),  as  a  function  of  X  for  all  E 
of  interest.  For  low  energies,  Eqs.  (92)  and  (93)  are  better 
converged  as  a  function  of  X  than  the  9  and  E  curves  indicate, 
due  to  the  convergence  of  gN^  itself  as  a  function  of  X.  The 
total  distribution  functions  for  both  bands  are  displaced  over¬ 
all  in  the  FxB  direction;  i.e.,  in  the  direction  of  the  Hall 
field.  This  was  ascertained  by  integrating  Eq.  (10b)  for  ij;^ 
over  the  hemisphere  in  the  direction  of  the  Hall  field. 
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Fig.  18.  Cubic  harmonic  decomposition  of  the  perturbed  part  6 
of  the  distribution  function,  Eq .  (74)  as  a  function 
of  hole  energy,  for  the  light  and  heavy  hole  bands 
of  silicon. 


V 


*«) 3  x-Ol 


E  ( 0.0045  eV) 


Fig.  19.  Cubic  harmonic  decomposition  of  the  perturbed  part  £ 
of  the  distribution  function,  Ea.  (77) ,  as  a  function 
of  hole  energy,  for  the  light  hole  band  of  silicon. 
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Pig.  20.  Cubic  harmonic  decomposition  of  the  perturbed  part  £ 
of  the  distribution  function,  Eq .  (77)  ,  as  a  function 
of  hole  energy,  for  the  heavy  hole  band  of  silicon. 
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Figs.  21-23  exhibit  functions  U^,  D^,  and  which  were 
defined  in  Eqs.  (79),  (80),  and  (81),  respectively.  The  plots 
for  UN,  Fig.  21,  display  marked  deviations  from  the  E  behavior 
expected  for  parabolic  bands.  The  orderina  U  >U  >U,,  indicates 
that  overall  the  H  band  will  dominate  the  total  conductivity 
mobility  although  its  partial  mobility  yc  is  lower  than  uc> 

The  spin-orbit  band  will  have  a  small  effect  on  the  mobility  at 
low  temperatures.  The  density  of  states  factors  DM,  Fig.  22, 
start  out  approximately  like  E  but  then  depart  considerably 
from  that  behavior.  The  ordering  DH>DL>Dg  is  expected  from  the 
band  curvatures.  A  most  complicated  behavior  for  WN  is  encountered 
for  silicon  in  Fig.  23.  The  E^  dependence  is  seen  mostly  at  the 
very  lowest  of  energies  examined.  Apparently,  the  approximate 
energy  dependences  of  various  factors  in  the  parabolic  limit, 
derived  in  Eqs.  (69-73)  and  (82-84),  break  down  completely  when 
compounded  to  produce  the  W  curves  for  silicon.  A  stronger  or 

weaker  energy  dependence  in  either  S,  G,  0 ,  H,  or  x,  from  that 

u 

expected  in  the  parabolic  limit,  would  produce  the  non-E  behavior 
of  WN.  This  is,  of  course,  expected  to  be  the  case  for  silicon 
with  its  nonparabolic  bands.  Fig.  23  indicates  that  the  light 
holes  in  the  L  band,  in  spite  of  their  lower  density  at  any  given 
temperature,  will  contribute  the  most  to  the  total  Hall  mobility. 

Figs.  24  and  25  display  the  calculated  partial  and  total 
conductivity  and  Hall  mobilities,  respectively.  As  expected,  the 
L  band  has  the  more  mobile  holes  and  dominates  the  overall  Hall 
mobility.  The  spin-orbit  band  makes  a  negligible  contribution  to 
either  calculated  mobility.  Fig.  26  presents  the  calculated  tem¬ 
perature  exponents  for  both  mobilities  and  the  r-factor  for 

silicon.  The  a  coefficient  starts  out  at  10K  with  a  vi.57  and 
c  c 

after  a  rapid  rise  saturates  in  value  at  about  2.10  for  T>90K. 

The  rise  of  from  1.55  at  10K  is  monotonic  in  the  temperature 
range  examined.  Therefore,  it  is  only  at  the  lowest  temperatures 
that  the  mobilities  come  close  to  the  parabolic  band  yvT  1-5 
behavior.  Also  shown  in  Fig.  26  are  temperature  exponents 
ac%l . 75+0 . 05  and  . 65+0 . 05  measured  by  us  in  the  20-40K  range. 
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Fig.  22.  Temperature  independent  part  of  the  integrand  in  the 
denominators  of  Eqs .  (53a  and  b)  for  u  and  uH, 
defined  by  Eq.  (80) . 
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Fig.  23.  Temperature  independent  part  of  the  integrand  in  the 
numerator  of  Eq.  (53b)  for  uH/  defined  by  Eq.  (81). 
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Fig.  24.  Calculated  total  and  partial  conductivity  mobilities, 
in  the  acoustic  phonon  limited  regime,  as  a  function 
of  temperature. 
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Fig.  25.  Calculated  total  and  partial  Hall  mobilities, 
in  the  acoustic  phonon  limited  regime,  as  a 
function  of  temperature. 
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Fig.  26.  Same  as  Fig.  10  for  silicon.  Experimental  data  for  a 
(open  circle)  and  aH  (filled  circle).  r-factor  data  from  ' 
sample  1202-H,  triangles  (Hall  bar  acceptors  6 . 57xl011cm~3 , 
donors  3 . 96x10 1 1^m~3)  and  sample  1300-V,  circles  (van  der  Pauw 
configuration,  acceptors  9 . 14xl011cm“3 ,  donors  3 . 33xl011cm-3 ) . 
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Figt  2G  shows  a  generally  close  agreement  between  experi¬ 
ment  and  theory  for  the  exponents.  Since  only  one  exponent  was 
fitted  to  the  experimental  data  it  is  difficult  to  say  whether 
the  exponents  have  an  appreciable  temperature  dependence  them¬ 
selves.  Measurements  encounter  difficulties  at  lower  temperatures, 
T<20K,  due  to  the  lack  of  signal  connected  with  a  low  thermal  hole 
population  in  high  purity  samples.  Using  lower  purity  samples  com¬ 
plicates  matters  due  to  the  presence  of  ionized  impurity  scatter¬ 
ing.  At  higher  temperatures,  T>40K,  the  problem  is  the  attain¬ 
ability  of  high  enough  magnetic  fields  to  saturate  the  Hall 
coefficient.  One  must  remember  that  this  saturation  is  more 
difficult  to  achieve  for  heavy  holes  due  to  their  lower  mobility. 

On  the  other  hand,  the  large  light  hole  mobilities  impose  more 
stringent  requirements  on  the  electric  fields  in  order  to  stay 
within  the  ohmic  regime. 

The  calculated  and  measured  r-factors  are  also  displayed  in 
Fig.  26.  The  two  sets  of  data  correspond  to  sample  1202-H  (Hall 
bar  configuration,  N.=6 . 57X101 1cm"3 ,  N=3. 96xl011cm~3)  and  sample 

^  ^  11-3 

1300-V  (van  der  Pauw  configuration,  N  =9. 14x10  cm  , 

11-3  A 

Nd=3. 33x10  cm  ).  Considering  the  spread  of  experimental  data 
the  agreement  between  experiment  and  theory  is  very  good  and  is 
on  the  order  of  4%. 

Fig.  27  shows  the  comparison  between  the  calculated  and 

measured  conductivity  mobility  (sample  1202-H) .  The  agreement  is 

seen  to  be  very  good.  There  is  a  hint  at  lower  temperatures  of 

the  influence  of  ionized  impurity  scattering  in  the  data  since  the 

first  two  points  fall  slightly  below  the  calculated  curve.  If 

this  conjecture  were  to  be  true  it  would  change  the  exponent  :<c 

in  Fig.  26  and  make  it  more  positive,  in  better  agreement  with 

the  calculated  a  . 

c 

The  low  field  limit  of  the  Hall  mobility  as  measured  by  us 
and  Elstner^66',  together  with  the  calculated  Hall  mobility  are 
shown  in  Fig.  28.  Again,  the  agreement  found  here  is  very  good. 

It  is  expected  that  at  about  100K  the  influence  of  the  optical 
phonon  scattering  will  be  felt.  The  theory  presented  in  this 
paper  can  be  extended  to  include  that  scattering  mechanism  as  well, 
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as  was  done  for  Ge  by  Lawaetz(4).  Note  that  these  two  phonon 
scattering  mechanisms  are  in  competition  with  one  another.  There¬ 
fore,  our  S  matrix,  Eq.  24,  will  consist  of  two  parts 


s  =  s  +  s 

ac  opt 


for  the  acoustic  and  optical  phonons, 
mobilities  we  need  to  invert  S;  e.g.. 


(92) 

Clearly,  to  solve  for  the 
Eqs .  (29)  and  (43), 


S'1  =  (S  +  S  .  ) 
ac  opt 


which  is  not  the  same  as 

„-l 


=  S  _1  +  S 
ac  opt 


-1 


-1 


(93) 


(94) 


This  last  generally  unacceptable  approximatior  leads  to  the  often 
used  expression 


u 


-1 


(95) 


l 

where  the  inverse  of  mobilities  for  each  scattering  mechanism  are 
added.  We  shall  summarize  results  of  this  calculation  and  measure¬ 
ments  in  the  next  section. 
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SECTION  VI 

SUMMARY  AND  CONCLUSIONS 


We  have  shown  that  the  experimental  Hall  factor  for  acoustic 
phonon  scattering  in  p-type  silicon  is  very  close  to  the  simple, 
textbook  value  of  3it/8.  There  is,  however,  a  slight  temperature 
dependence  that  is  not  predicted  by  theories  that  neglect  the 
nonparabolic  nature  of  the  valence  bands.  This  temperature 
dependent  Hall  factor  has  been  shown  to  have  a  noticeable  effect 
on  the  results  of  Hall  effect  analysis,  affecting  primarily  the 
acceptor  concentration.  We  used  the  experimentally  determined 
r's  to  calculate  the  true  conductivity  mobility  for  acoustic 
phonon  scattering  as  a  function  of  temperature.  These  values 
were  in  turn  compared  with  our  calculations  where  agreement  was 
found  to  be  excellent.  The  experiment  yields  a  temperature 
dependence  of  T  1,75  for  u  and  T-1'65  for  u„(H-*-0). 

C  ri 

The  conductivity  and  Hall  mobilities  for  silicon  and 
germanium  have  been  calculated  using  the  deformation  potential 
theory  of  Tiersten(3)  and  the  formalism  developed  by  Lawaetz<4). 

We  have  solved  the  Boltzmann  equation  without  recourse  to  the 
relaxation  time  approximation.  The  agreement  between  experimental 
and  theoretical  r-factors  is  on  the  order  of  4%  for  silicon  and 
12%  for  germanium.  The  corresponding  mobilities  are  generally 
in  better  agreement  with  the  experiments  than  are  the  r-factors. 

The  results  of  our  calculation  of  the  conductivity  and 
Hall  mobilities  for  Si  and  Ge,  together  with  the  pioneering 
works  of  Tiersten  and  Lawaetz,  indicate  that  phonon  limited 
transport  in  nonpolar  semiconductors  can  be  quantitatively 
modeled  within  the  framework  of  the  first-principles  deformation 
potential  theory.  It  appears  that  the  main  limitation  to  the 
quantitative  accuracy  of  the  method  is  the  precision  with  which 
the  deformation  potentials  are  known,  other  input  parameters 
having  already  been  determined  with  sufficient  precision. 

Elimination  of  a  few  approximations  in  our  theoretical 
treatment  can  make  the  theoretical  results  more  exact  but  that 


75 


might  be  beyond  the  resolution  of  present  day  experiments. 

These  approximations,  taken  together  with  the  inherent  limita¬ 
tions  of  the  deformation  potential  concept  itself,  can  be  esti¬ 
mated  to  produce  a  +10%  fundamental  accuracy  of  the  method 
employed  here(3'4*.  These  approximations  will  produce  worse 
agreement  for  progressively  higher  transport  coefficients;  e.g., 
longitudinal  and  transverse  magnetoresistance,  etc. 

The  lowest  temperature  results  would  be  improved  by 
retaining  the  phonon  energy  in  the  delta  function  argument, 

Eq.  (3).  Another  improvement  would  result  from  the  use  of  the 
dynamical  matrix  solutions  for  phonon  energies  and  polarization 
vectors.  By  using  a  degenerate  perturbation  theory,  one  also 
should  eliminate  the  ambiguity  in  defining  the  polarization 
vectors  in  the  q-+0  limit  for  forward  intraband  hole  scatter¬ 
ing^3'4^.  Finally,  as  Lawaetz  has  shown^G7^,  the  deformation 
potential  scattering  theory  is  only  the  first  order  long  wave¬ 
length  approximation  to  the  full  carrier-phonon  interaction  term. 
This  approximation  he  shows  is  good  to  at  least  5%,  which  pro¬ 
vides  the  limit  on  the  accuracy  of  the  method  itself. 
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APPENDIX  A 

SYMMETRY  ANALYSIS  FOR  CONDUCTIVITY  TENSORS 


The  most  general  form  for  the  current  density 


is 


(68) 


J  —  0  F  +0  ^  F  B-.+0  ..  F  B ,  B  + 

u  py  v  yvA  v  A  yvAn  v  a  n 


(Al) 


to  terms  linear  in  the  electric  field  F,  where  the  sums  are  to 
be  performed  over  repeated  indices  which  stand  for  x,  y,  and  z. 
The  first  term  in  the  sum  can  be  written  in  terms  of  the  con- 
ductivity  tensor  a  as 

J  =  "o’  •  F  ,  (A2 ) 

where  both  J  and  F  transform  like  vectors 


J  =  J' 


F  =  F' 


(A3) 


when  operated  on  by  one  of  the  48  rotations  R  of  the  cubic  point 
group  0,.  In  the  rotated  coordinate,  system  the  conductivity  tensor 
a  becomes  o',  where 


so  that  from  (A2), 


J'  =  o'  • 

(A3),  and  (A4) 

*-*■.  *-*■ 
o'  =  R  • 


F' 

(69) 


-1 


R 


( A4 ) 


(A5 ) 


for  all  ReO^.  Since  a  is  a  macroscopic  property  of  the  material, 
it  is  invariant  under  all  rotations,  so  that 

*o'  =  V  .  ( A6 ) 


Applying  a  few  rotation  matrices  to  the  tensor  o  in  Eq.  (A5) 
quickly  establishes  that 


and 


a  =0  for  y  t*  v 

yv 


a 

o 


(A7a ) 


so  that  the  second  rank  conductivity  tensor  is  diagonal . 

(69) 

A  similar  analysis  for  the  third  rank  tensor  yields 

Jy  ~  ayvXFvRX  ' 

RJy  =  °ivA(BPv)(RBX)  • 


(A8 ) 
(A9 ) 
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Vectors  J  and  F  transform  according  to  the  rows  of  the 
representation  whose  matrix  representation  is  given  by  the 
rotation  matrices  R  themselves.  The  magnetic  field  B  trans¬ 
forms  like  a  pseudovector  under  rotations  according  to  the  rows 
of  the  T'5  representation.  The  matrix  representation  for  F^ 
is  given  by  the  rotation  matrices  R  for  the  proper  rotations  and 
(-IR)  for  improper  rotations,  where  I  is  the  inversion  matrix. 
Again,  the  a'  tensor  is  invariant  under  rotations  to  that 


a '  ,  =  o  , 
y  v  A  yvA 

which  leads  via  Eqs.  (A8)  and  (A9)  to  the  condition 


(A10) 


yv  A 


R  , R  ,F™  ,(R)0  ,  , 

y  y  vv  A  A  y  v  A 


y '  v  '  A  ' 


where  Fm(R)  are  the  matrices  for  the  representation.  Per¬ 

forming  the  indicated  operations  in  Eq.  (All) ,  with  a  few  rota¬ 
tion  operators,  shows  that 


ayvA  €yvAa 

where  €  ,  is  the  Levi-Cevita  symbol. 

V  A 

to  first  order  in  $  and  B 


(All) 

Altogether,  we  obtain 


Clearly,  a 

-> 

J 


so  that 


J  =  a  F+a(FxB) 
c 

epuc,  and 

ep<V>  =  ep  (<V^g>F+<V^>FB+  ...)  , 


y 


c 


•  F 


9 


OL 


FxB 


|  FxB | 2 


(A12 ) 


(A13) 

(A14 ) 

( A15 ) 


The  coefficient  a  is  related  to  the  Hall  mobility,  yH. 
mobility  is  defined  as 


UH 


(eVh 

BzFx 


The  Hall 


(A16 ) 
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where  Fx  is  the  electric  field  in  direction  of  the  applied 

voltage,  (F  )H  the  field  in  direction  of  the  Hall  voltage,  and 

B  is  the  applied  magnetic  field.  More  generally,  the  Hall 
z 

electric  field  is  given  by 

Fh  =  uh(FxB)  ,  (A17) 

so  that 

J  =  epuc[F+FHJ 

=  ep  [V‘cF+ycU{J  (?xBJ  ]  .  (A18) 

By  comparison  with  (A12) 

2 

a  =  ucUH  =  yc(ruc)  =  ruc  .  (A19) 

It  should  be  said  that  the  above  proof  implicitly  uses  the 
Neumann  theorem  which  states  that  macroscopic  physical  properties 
of  a  crystal  are  invariant  under  point  group  operations  of  the 
crystal's  soace  group.  In  case  of  the  diamond  structure,  fraction¬ 
al  translations  in  a  unit  cell  are  too  small  on  the  macroscopic 
scale  to  influence  the  results  of  experimental  measurements. 
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APPENDIX  B 


COUPLING  COEFFICIENTS  FOR  THE 
MAGNETIC  FIELD  TERM,  EQ .  64 


An  alternative  technique  for  dealing  with  the  operations 
indicated  by  Eq.  (64)  are  provided  in  this  Appendix.  It  is 
necessary  to  work  with  the  expression 

hN(E,k)K^(k)  =  t%EN(k)x7j^]EKjM(k)  .  (B1  ] 

We  recast  the  energy  gradient  into  the  form  of  Eq.  (62b)  and 
obtain 


2 Y"*2  €  (x  IT.+x3 
N  qpr  Ip  4  p 


h",E'*)KN>  ■  \Pr 

[(Xp'I4+!tpII5+xpn6)  SIT:] 


*N>  * 


*{vlk)  , 


where  xg  are  direction  cosines,  xq=kq/k .  Typical  terms  in 
Eq.  (B2)  are 


= 

q  x 

_€qpr 

XP 

3Kr] 

K\ 

9 

(B3a) 

,(3)„6u  _ 

■q  KX  = 

,€qpr 

x3 

P 

£] 

Kx 

9 

(B3b) 

(5)„Sy  = 
q  KX  ~ 

L€qpr 

5 

X 

P 

9 

3kr- 

Kx 

* 

(B3c) 

The  operations  in  Eqs.  (B3)  can  be  expressed  as  a  sum  of  cubic 
harmonics  with  coupling  coefficients  defined  by 

■“W-S  c^Vx,k»'  • 

We  have  worked  out  a  number  of  tables  for  the  coupling 
coefficients  C  defined  by  Eq.  (B4) .  Some  of  the  tables  are 
provided  below  in  the  VdLB  notation ( 46) ; 


kki1'* 


q  ip 


o 


We  note  that  these  matrices  are  traceless  and  that  the  matrices 
with  5's  inside  are  antisymmetric,  while  those  with  e'  are  symmetric. 
Although  we  have  a  number  of  other  matrices  for  these  tables,  our 
tables  for  are  not  complete  enough  to  warrant  their 

*  "i 

reporting.  We  would  encourage  other  workers  to  complete  these 
tables.  They  have  the  potential  of  simplifying  the  computational 
effort  in  an  elegant  manner.  This  will  be  especially  true  in  cal¬ 
culations  of  higher  magnetoresistance  coefficients. 
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